We use recent 36 observational Hubble data (OHD) in the redshift range $0.07\leq z\leq 2.36$, latest \textgravedbl joint light curves\textacutedbl (JLA) sample, comprised of 740 type Ia supernovae (SNIa) in the redshift range $0.01\leq z \leq 1.30$, and their joint combination datasets to constrain anisotropic Bianchi type I (BI) dark energy (DE) model. To estimate model parameters, we apply Hamiltonian Monte Carlo technique. We also compute the covariance matrix for BI dark energy model by considering different datasets to compare the correlation between model parameters. To check the acceptability of our fittings, all results are compared with those obtained from 9 year WMAP as well as Planck (2015) collaboration. Our estimations show that at 68\% confidence level the dark energy equation of state (EOS) parameter for OHD or JLA datasets alone varies between quintessence and phantom regions whereas for OHD+JLA dataset this parameter only varies in phantom region. It is also found that the current cosmic anisotropy is of order $\sim10^{-3}$ which imply that the OHD and JLA datasets do not put tight constraint on this parameter. Therefore, to constraint anisotropy parameter, it is necessary to use high redshif dataset namely cosmic microwave background (CMB). Moreover, from the calculation of $p$-value associated with $\chi^{2}$ statistic we observed that non of the $\omega \mbox{BI}$ and flat $\omega\mbox{CDM}$ models rule out by OHD or JLA datasets. The deceleration parameter is obtained as $q=-0.46^{+0.89 +0.36}_{-0.41 -0.37}$, $q=-0.619^{+0.12 +0.20}_{-0.095 -0.24}$, and $q=-0.52^{+0.080 +0.014}_{-0.046 -0.15}$ for OHD, SNIa, and OHD+SNIa data respectively.

Using the mathematical definitions of deceleration and jerk parameters we obtain a general differential equation for squared Hubble parameter. For a constant jerk, this differential equation leads to an exact function for Hubble parameter. By the aid of this exact Hubble function we can exactly reconstruct any other cosmographic parameters. We also obtained a general function for transition redshift as well as spacetime curvature. Our derived functions clearly impose a lower limit on the jerk parameter which is $j_{min}\geq-0.125$. Moreover, we found that the jerk parameter indicates the geometry of the spacetime i.e any deviation from $j=1$ imply to a non-flat spacetime. In other word $j\neq 1$ reefers to a dynamical, time varying, dark energy. From obtained Hubble function we recover the analogue of $\Lambda$CDM model. To constrain cosmographic parameters as well as transition redshift and spacetime curvature of the recovered $\Lambda$CDM model, we used Metropolis-Hasting algorithm to perform Monte Carlo Markov Chain analysis by using observational Hubble data obtained from cosmic chronometric (CC) technique, BAO data, Pantheon compilation of Supernovae type Ia, and their joint combination. The only free parameters are $H$, $A(\Omega_{m})$ and $j$. From joint analysis we obtained $H_{0}=69.9\pm 1.7$, $A(\sim\Omega_{0m})=0.279^{+0.013}_{-0.017}$, $B(\sim\Omega_{0X})=0.721^{+0.017}_{-0.013}$, $j_{0}=1.038^{+0.061}_{-0.023}$ and $z_{t}=0.706^{+0.031}_{-0.034}$.

We use Gaussian processes in combination with MCMC method to place constraints on cosmological parameters of three dark energy models including flat and curved FRW and Bianchi type I spacetimes. To do so, we use recently compiled 36 measurements of the Hubble parameter $H(z)$ in the redshifts intermediate $0.07\leqslant z \leqslant 2.36$. Moreover, we use these models to estimate the redshift of the deceleration-acceleration transition. We consider two Gaussian priors for current value of the Hubble constant i.e $H_{0}=73\pm1.74 (68\pm 2.8)$ km/s/Mpc to investigate the effect of the assumed $H_{0}$ on our parameters estimations. For statistical analysis we use NUTS sampler which is an extension of Hamiltonian Monte Carlo algorithm to generate MCMC chains for parameters of dark energy models. To compare the considered cosmologies, we perform Akaike information criterion (AIC) and Bayes factor ($\Psi$). In general, when we compared our results with 9 years WMAP as well as Planck 2015 Collaboration, we found that Bianchi type I model is slightly fits better to the observational Hubble data with respect to the non-flat FRW model.