From 082d58a7ca7f6b833b13a50746f50907a2dd473c Mon Sep 17 00:00:00 2001
From: "github-actions[bot]"
 <41898282+github-actions[bot]@users.noreply.github.com>
Date: Thu, 3 Jun 2021 02:42:48 +0000
Subject: [PATCH 1/6] CompatHelper: bump compat for "Turing" to "0.16"

---
 Project.toml | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/Project.toml b/Project.toml
index 4d734278..246b7943 100644
--- a/Project.toml
+++ b/Project.toml
@@ -53,7 +53,7 @@ Reexport = "0.2, 1.0"
 Requires = "0.5, 1.0"
 StructArrays = "0.4, 0.5"
 TransformVariables = "0.3, 0.4"
-Turing = "0.12, 0.13, 0.14, 0.15"
+Turing = "0.12, 0.13, 0.14, 0.15, 0.16"
 julia = "1.3"
 
 [extras]

From cf92af6050c289bfec652dcaa6cf7fa6185c8cb7 Mon Sep 17 00:00:00 2001
From: Vaibhav Dixit <vaibhavyashdixit@gmail.com>
Date: Fri, 4 Jun 2021 16:33:53 +0530
Subject: [PATCH 2/6] Update `stan_inference` to new ode syntax of stan

---
 src/stan_inference.jl | 27 ++++++++++++++-------------
 1 file changed, 14 insertions(+), 13 deletions(-)

diff --git a/src/stan_inference.jl b/src/stan_inference.jl
index 126c2a0f..f17692cc 100644
--- a/src/stan_inference.jl
+++ b/src/stan_inference.jl
@@ -10,13 +10,13 @@ end
 
 function generate_priors(n,priors)
   priors_string = ""
-  if priors==nothing
+  if priors===nothing
     for i in 1:n
-      priors_string = string(priors_string,"theta[$i] ~ normal(0, 1)", " ; ")
+      priors_string = string(priors_string,"theta_$i ~ normal(0, 1)", " ; ")
     end
   else
     for i in 1:n
-      priors_string = string(priors_string,"theta[$i] ~",stan_string(priors[i]),";")
+      priors_string = string(priors_string,"theta_$i ~",stan_string(priors[i]),";")
     end
   end
   priors_string
@@ -66,11 +66,11 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
   
   if stanmodel === nothing
     if alg ==:adams
-      algorithm = "integrate_ode_adams"
+      algorithm = "ode_adams_tol"
     elseif alg ==:rk45
-      algorithm = "integrate_ode_rk45"
+      algorithm = "ode_rk45_tol"
     elseif alg == :bdf
-      algorithm = "integrate_ode_bdf"
+      algorithm = "ode_bdf_tol"
     else
       error("The choices for alg are :adams, :rk45, or :bdf")
     end
@@ -78,9 +78,11 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
     tuple_hyper_params = ""
     setup_params = ""
     thetas = ""
+    theta_names = ""
     theta_string = generate_theta(length_of_parameter,priors)
     for i in 1:length_of_parameter
-      thetas = string(thetas,"theta[$i] = theta$i",";")
+      thetas = string(thetas,"real theta_$i",";")
+      theta_names = string(theta_names,"theta_$i",",")
     end
     for i in 1:length_of_params
       if isa(vars[i],StanODEData)
@@ -103,12 +105,12 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
           u0 = u0*string(u_)
         end
         u0 = u0*"}"
-        integral_string = "u_hat = $algorithm(sho, append_array(theta[1:$nu],$u0), t0, ts, theta[$(nu+1):$length_of_parameter], x_r, x_i, $reltol, $abstol, $maxiter);"
+        integral_string = "u_hat = $algorithm(sho, append_array(theta[1:$nu],$u0), t0, ts, $reltol, $abstol, $maxiter, theta[$(nu+1):$length_of_parameter]);"
       else 
-        integral_string = "u_hat = $algorithm(sho, theta[1:$nu], t0, ts, theta[$(nu+1):$length_of_parameter], x_r, x_i, $reltol, $abstol, $maxiter);"
+        integral_string = "u_hat = $algorithm(sho, theta[1:$nu], t0, ts, $reltol, $abstol, $maxiter, theta[$(nu+1):$length_of_parameter]);"
       end
     else
-      integral_string = "u_hat = $algorithm(sho, u0, t0, ts, theta, x_r, x_i, $reltol, $abstol, $maxiter);"
+      integral_string = "u_hat = $algorithm(sho, u0, t0, ts, $reltol, $abstol, $maxiter, $(rstrip(theta_names,',')));"
     end
     binsearch_string = """
       int bin_search(real x, int min_val, int max_val){
@@ -141,7 +143,7 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
       $diffeq_string
     }
     data {
-      real u0[$length_of_y];
+      vector[$length_of_y] u0;
       int<lower=1> T;
       real internal_var___u[T,$(length(save_idxs))];
       real t0;
@@ -156,11 +158,10 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
       $theta_string
     }
     transformed parameters{
-      real theta[$length_of_parameter];
       $thetas
     }
     model{
-      real u_hat[T,$length_of_y];
+      vector[$length_of_y] u_hat[T];
       $hyper_params
       $priors_string
       $integral_string

From fcaf3a3ee0c039ce6724a08bfbeb783b89d26296 Mon Sep 17 00:00:00 2001
From: Vaibhav Dixit <vaibhavyashdixit@gmail.com>
Date: Sat, 10 Jul 2021 21:02:54 +0530
Subject: [PATCH 3/6] Update stan_inference differential equation string

---
 Project.toml          |  2 +-
 src/DiffEqBayes.jl    | 14 ++++----------
 src/stan_inference.jl | 30 +++++++++++-------------------
 3 files changed, 16 insertions(+), 30 deletions(-)

diff --git a/Project.toml b/Project.toml
index 246b7943..59e04c04 100644
--- a/Project.toml
+++ b/Project.toml
@@ -26,6 +26,7 @@ Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 Requires = "ae029012-a4dd-5104-9daa-d747884805df"
+StanSample = "c1514b29-d3a0-5178-b312-660c88baa699"
 StructArrays = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
 TransformVariables = "84d833dd-6860-57f9-a1a7-6da5db126cff"
 Turing = "fce5fe82-541a-59a6-adf8-730c64b5f9a0"
@@ -46,7 +47,6 @@ Missings = "0.4, 1.0"
 ModelingToolkit = "5.6"
 Optim = "0.19, 0.20, 0.21, 0.22, 1.0"
 PDMats = "0.9, 0.10, 0.11"
-ParameterizedFunctions = "5"
 Parameters = "0.12"
 RecursiveArrayTools = "1,2"
 Reexport = "0.2, 1.0"
diff --git a/src/DiffEqBayes.jl b/src/DiffEqBayes.jl
index 6cb85ecb..7c662b67 100644
--- a/src/DiffEqBayes.jl
+++ b/src/DiffEqBayes.jl
@@ -7,22 +7,17 @@ using DocStringExtensions
 using DiffEqBase, Distributions, Turing, MacroTools
 using RecursiveArrayTools, ModelingToolkit
 using Parameters, Distributions, Optim, Requires
-using Distances, ApproxBayes, DocStringExtensions, Random
+using Distances, ApproxBayes, DocStringExtensions, Random, StanSample
 
 STANDARD_PROB_GENERATOR(prob,p) = remake(prob;u0=eltype(p).(prob.u0),p=p)
 STANDARD_PROB_GENERATOR(prob::EnsembleProblem,p) = EnsembleProblem(remake(prob.prob;u0=eltype(p).(prob.prob.u0),p=p))
 
 include("turing_inference.jl")
 include("abc_inference.jl")
+include("stan_string.jl")
+include("stan_inference.jl")
 
 function __init__()
-    @require CmdStan="593b3428-ca2f-500c-ae53-031589ec8ddd" begin
-        using .CmdStan
-        include("stan_inference.jl")
-        include("stan_string.jl")
-        export stan_inference, stan_string
-    end
-
     @require DynamicHMC="bbc10e6e-7c05-544b-b16e-64fede858acb" begin
         using .DynamicHMC, TransformVariables, LogDensityProblems
         include("dynamichmc_inference.jl")
@@ -30,6 +25,5 @@ function __init__()
     end
 end
 
-export turing_inference, abc_inference
-
+export turing_inference, stan_inference ,abc_inference
 end # module
diff --git a/src/stan_inference.jl b/src/stan_inference.jl
index f17692cc..4a7f1b61 100644
--- a/src/stan_inference.jl
+++ b/src/stan_inference.jl
@@ -16,7 +16,7 @@ function generate_priors(n,priors)
     end
   else
     for i in 1:n
-      priors_string = string(priors_string,"theta_$i ~",stan_string(priors[i]),";")
+      priors_string = string(priors_string,"theta_$i ~ ",stan_string(priors[i]),";")
     end
   end
   priors_string
@@ -34,13 +34,13 @@ function generate_theta(n,priors)
       lower_bound = string("lower=",minimum(priors[i]))
     end
     if lower_bound != "" && upper_bound != ""
-      theta = string(theta,"real","<$lower_bound",",","$upper_bound>"," theta$i",";")
+      theta = string(theta,"real","<$lower_bound",",","$upper_bound>"," theta_$i",";")
     elseif lower_bound != ""
-      theta = string(theta,"real","<$lower_bound",">"," theta$i",";")
+      theta = string(theta,"real","<$lower_bound",">"," theta_$i",";")
     elseif upper_bound != ""
-      theta = string(theta,"real","<","$upper_bound>"," theta$i",";")
+      theta = string(theta,"real","<","$upper_bound>"," theta_$i",";")
     else
-      theta = string(theta,"real"," theta$i",";")
+      theta = string(theta,"real"," theta_$i",";")
     end
   end
   return theta
@@ -63,7 +63,7 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
   else
     length_of_parameter = length(prob.p) + sample_u0 * length(save_idxs)
   end
-  
+
   if stanmodel === nothing
     if alg ==:adams
       algorithm = "ode_adams_tol"
@@ -106,7 +106,7 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
         end
         u0 = u0*"}"
         integral_string = "u_hat = $algorithm(sho, append_array(theta[1:$nu],$u0), t0, ts, $reltol, $abstol, $maxiter, theta[$(nu+1):$length_of_parameter]);"
-      else 
+      else
         integral_string = "u_hat = $algorithm(sho, theta[1:$nu], t0, ts, $reltol, $abstol, $maxiter, theta[$(nu+1):$length_of_parameter]);"
       end
     else
@@ -122,8 +122,8 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
                 out = mid_pt;
                 range = 0;
             } else {
-                range = (range + 1) / 2; 
-                mid_pt = x > mid_pt ? mid_pt + range: mid_pt - range; 
+                range = (range + 1) / 2;
+                mid_pt = x > mid_pt ? mid_pt + range: mid_pt - range;
             }
         }
         return out;
@@ -149,17 +149,10 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
       real t0;
       real ts[T];
     }
-    transformed data {
-      real x_r[0];
-      int x_i[0];
-    }
     parameters {
       $setup_params
       $theta_string
     }
-    transformed parameters{
-      $thetas
-    }
     model{
       vector[$length_of_y] u_hat[T];
       $hyper_params
@@ -170,9 +163,8 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
         }
     }
     "
-    stanmodel = CmdStan.Stanmodel(num_samples=num_samples, num_warmup=num_warmup, name="parameter_estimation_model", model=parameter_estimation_model, nchains=nchains, printsummary = printsummary)
+    stanmodel = StanSample.SampleModel("parameter_estimation_model", parameter_estimation_model)
   end
   parameter_estimation_data = Dict("u0"=>prob.u0, "T" => length(t), "internal_var___u" => view(data, :, 1:length(t))', "t0" => prob.tspan[1], "ts" => t)
-  return_code, chains, cnames = CmdStan.stan(stanmodel, [parameter_estimation_data])
-  return StanModel(stanmodel, return_code, chains, cnames)
+  stan_sample(stanmodel; data = parameter_estimation_data)
 end

From 06f9e33dcc641cc4037efa90642ec3c82630742e Mon Sep 17 00:00:00 2001
From: Vaibhav Dixit <vaibhavyashdixit@gmail.com>
Date: Sun, 11 Jul 2021 12:51:22 +0530
Subject: [PATCH 4/6] Update tests and results are returned in a struct

---
 Project.toml          |  2 --
 src/DiffEqBayes.jl    |  4 ++--
 src/stan_inference.jl | 31 ++++++++++++++++++++-----------
 test/runtests.jl      |  4 +---
 test/stan.jl          | 32 +++++++++++++-------------------
 5 files changed, 36 insertions(+), 37 deletions(-)

diff --git a/Project.toml b/Project.toml
index 59e04c04..72083ba7 100644
--- a/Project.toml
+++ b/Project.toml
@@ -4,7 +4,6 @@ authors = ["Vaibhavdixit02 <vaibhavyashdixit@gmail.com>"]
 version = "2.24.0"
 
 [deps]
-ApproxBayes = "f5f396d3-230c-5e07-80e6-9fadf06146cc"
 DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
 DiffResults = "163ba53b-c6d8-5494-b064-1a9d43ac40c5"
 Distances = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7"
@@ -32,7 +31,6 @@ TransformVariables = "84d833dd-6860-57f9-a1a7-6da5db126cff"
 Turing = "fce5fe82-541a-59a6-adf8-730c64b5f9a0"
 
 [compat]
-ApproxBayes = "0.3"
 DiffEqBase = "6.36"
 DiffResults = "0.0.4, 1.0"
 Distances = "0.8, 0.9, 0.10"
diff --git a/src/DiffEqBayes.jl b/src/DiffEqBayes.jl
index 7c662b67..3df2043a 100644
--- a/src/DiffEqBayes.jl
+++ b/src/DiffEqBayes.jl
@@ -7,13 +7,13 @@ using DocStringExtensions
 using DiffEqBase, Distributions, Turing, MacroTools
 using RecursiveArrayTools, ModelingToolkit
 using Parameters, Distributions, Optim, Requires
-using Distances, ApproxBayes, DocStringExtensions, Random, StanSample
+using Distances, DocStringExtensions, Random, StanSample
 
 STANDARD_PROB_GENERATOR(prob,p) = remake(prob;u0=eltype(p).(prob.u0),p=p)
 STANDARD_PROB_GENERATOR(prob::EnsembleProblem,p) = EnsembleProblem(remake(prob.prob;u0=eltype(p).(prob.prob.u0),p=p))
 
 include("turing_inference.jl")
-include("abc_inference.jl")
+# include("abc_inference.jl")
 include("stan_string.jl")
 include("stan_inference.jl")
 
diff --git a/src/stan_inference.jl b/src/stan_inference.jl
index 4a7f1b61..c41903dc 100644
--- a/src/stan_inference.jl
+++ b/src/stan_inference.jl
@@ -1,8 +1,11 @@
-struct StanModel{M,R,C,N}
+struct StanResult{M,R,C}
   model::M
   return_code::R
   chains::C
-  cnames::N
+end
+
+function Base.show(io::IO, mime::MIME"text/plain", res::StanResult)
+  show(io, mime, res.chains)
 end
 
 struct StanODEData
@@ -50,9 +53,10 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
                             stanmodel = nothing;alg=:rk45,
                             num_samples=1000, num_warmup=1000, reltol=1e-3,
                             abstol=1e-6, maxiter=Int(1e5),likelihood=Normal,
-                            vars=(StanODEData(),InverseGamma(3,3)),nchains=1,
-                            sample_u0 = false, save_idxs = nothing, diffeq_string = nothing, printsummary = true)
-  
+                            vars=(StanODEData(),InverseGamma(3,3)),nchains=[1],
+                            sample_u0 = false, save_idxs = nothing, diffeq_string = nothing, 
+                            printsummary = true, output_format = :mcmcchains)
+
   save_idxs !== nothing && length(save_idxs) == 1 ? save_idxs = save_idxs[1] : save_idxs = save_idxs
   length_of_y = length(prob.u0)
   save_idxs = something(save_idxs, 1:length_of_y)
@@ -99,15 +103,15 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
     stan_likelihood = stan_string(likelihood)
     if sample_u0
       nu = length(save_idxs)
+      dv_names_ind = findfirst("$nu", theta_names)[1]
       if nu < length(prob.u0)
-        u0 = "{"
+        u0 = ""
         for u_ in prob.u0[nu+1:length(prob.u0)]
           u0 = u0*string(u_)
         end
-        u0 = u0*"}"
-        integral_string = "u_hat = $algorithm(sho, append_array(theta[1:$nu],$u0), t0, ts, $reltol, $abstol, $maxiter, theta[$(nu+1):$length_of_parameter]);"
+        integral_string = "u_hat = $algorithm(sho, [$(theta_names[1:dv_names_ind]),$u0]', t0, ts, $reltol, $abstol, $maxiter, $(rstrip(theta_names[dv_names_ind+2:end],',')));"
       else
-        integral_string = "u_hat = $algorithm(sho, theta[1:$nu], t0, ts, $reltol, $abstol, $maxiter, theta[$(nu+1):$length_of_parameter]);"
+        integral_string = "u_hat = $algorithm(sho, [$(theta_names[1:dv_names_ind])]', t0, ts, $reltol, $abstol, $maxiter, $(rstrip(theta_names[dv_names_ind+2:end],',')));"
       end
     else
       integral_string = "u_hat = $algorithm(sho, u0, t0, ts, $reltol, $abstol, $maxiter, $(rstrip(theta_names,',')));"
@@ -163,8 +167,13 @@ function stan_inference(prob::DiffEqBase.DEProblem,t,data,priors = nothing,
         }
     }
     "
-    stanmodel = StanSample.SampleModel("parameter_estimation_model", parameter_estimation_model)
+    stanmodel = StanSample.SampleModel("parameter_estimation_model", parameter_estimation_model, nchains; printsummary = printsummary, method = StanSample.Sample(;num_samples = num_samples, num_warmup = num_warmup))
   end
   parameter_estimation_data = Dict("u0"=>prob.u0, "T" => length(t), "internal_var___u" => view(data, :, 1:length(t))', "t0" => prob.tspan[1], "ts" => t)
-  stan_sample(stanmodel; data = parameter_estimation_data)
+  rc = stan_sample(stanmodel; data = parameter_estimation_data)
+  if success(rc)
+    return StanResult(stanmodel, rc, read_samples(stanmodel; output_format=output_format))
+  else
+     rc.err
+  end
 end
diff --git a/test/runtests.jl b/test/runtests.jl
index 4fb0e52f..2ac9295f 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -8,12 +8,10 @@ const GROUP = get(ENV, "GROUP", "All")
 if GROUP == "All" || GROUP == "Core"
     @time @safetestset "DynamicHMC" begin include("dynamicHMC.jl") end
     @time @safetestset "Turing" begin include("turing.jl") end
-    @time @safetestset "ABC" begin include("abc.jl") end
+    # @time @safetestset "ABC" begin include("abc.jl") end
 end
 
 if GROUP == "Stan" || GROUP == "All"
-    using Pkg
-    Pkg.add("CmdStan")
     @time @safetestset "Stan_String" begin include("stan_string.jl") end
     @time @safetestset "Stan" begin include("stan.jl") end
 end
diff --git a/test/stan.jl b/test/stan.jl
index b7e4ea77..280d6151 100644
--- a/test/stan.jl
+++ b/test/stan.jl
@@ -1,4 +1,4 @@
-using CmdStan, DiffEqBayes, OrdinaryDiffEq, ParameterizedFunctions,
+using DiffEqBayes, OrdinaryDiffEq, ParameterizedFunctions,
       RecursiveArrayTools, Distributions, Test
 
 println("One parameter case")
@@ -19,24 +19,21 @@ priors = [truncated(Normal(1.5,0.1),1.0,1.8)]
 bayesian_result = stan_inference(prob1,t,data,priors;num_samples=300,
                                  num_warmup=500,likelihood=Normal)
 
-sdf  = CmdStan.read_summary(bayesian_result.model)
-@test sdf[sdf.parameters .== :theta1, :mean][1] ≈ 1.5 atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_1)[1]) ≈ 1.5 atol=3e-1
 
 # Test norecompile
 bayesian_result2 = stan_inference(prob1,t,data,priors,bayesian_result.model;
                                   num_samples=300,num_warmup=500,likelihood=Normal)
 
-sdf  = CmdStan.read_summary(bayesian_result.model)
-@test sdf[sdf.parameters .== :theta1, :mean][1] ≈ 1.5 atol=3e-1
+@test mean(get(bayesian_result2.chains,:theta_1)[1]) ≈ 1.5 atol=3e-1
 
 priors = [truncated(Normal(1.,0.01),0.5,2.0),truncated(Normal(1.,0.01),0.5,2.0),truncated(Normal(1.5,0.01),1.0,2.0)]
 bayesian_result = stan_inference(prob1,t,data,priors;num_samples=300,
                                  num_warmup=500,likelihood=Normal,sample_u0=true)
 
-sdf  = CmdStan.read_summary(bayesian_result.model)
-@test sdf[sdf.parameters .== :theta1, :mean][1] ≈ 1. atol=3e-1
-@test sdf[sdf.parameters .== :theta2, :mean][1] ≈ 1. atol=3e-1
-@test sdf[sdf.parameters .== :theta3, :mean][1] ≈ 1.5 atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_1)[1]) ≈ 1. atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_2)[1]) ≈ 1. atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_3)[1]) ≈ 1.5 atol=3e-1
 
 sol = solve(prob1,Tsit5(),save_idxs=[1])
 randomized = VectorOfArray([(sol(t[i]) + .01 * randn(1)) for i in 1:length(t)])
@@ -45,17 +42,15 @@ priors = [truncated(Normal(1.5,0.1),0.5,2)]
 bayesian_result = stan_inference(prob1,t,data,priors;num_samples=300,
                                  num_warmup=500,likelihood=Normal,save_idxs=[1])
 
-sdf  = CmdStan.read_summary(bayesian_result.model)
-@test sdf[sdf.parameters .== :theta1, :mean][1] ≈ 1.5 atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_1)[1]) ≈ 1.5 atol=3e-1
 
 
 priors = [truncated(Normal(1.,0.01),0.5,2),truncated(Normal(1.5,0.01),0.5,2)]
 bayesian_result = stan_inference(prob1,t,data,priors;num_samples=300,
                                  num_warmup=500,likelihood=Normal,save_idxs=[1],sample_u0=true)
 
-sdf  = CmdStan.read_summary(bayesian_result.model)
-@test sdf[sdf.parameters .== :theta1, :mean][1] ≈ 1. atol=3e-1
-@test sdf[sdf.parameters .== :theta2, :mean][1] ≈ 1.5 atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_1)[1]) ≈ 1. atol=3e-1
+@test mean(get(bayesian_result.chains,:theta_2)[1]) ≈ 1.5 atol=3e-1
 
 println("Four parameter case")
 f1 = @ode_def begin
@@ -74,8 +69,7 @@ priors = [truncated(Normal(1.5,0.01),0.5,2),truncated(Normal(1.0,0.01),0.5,1.5),
           truncated(Normal(3.0,0.01),0.5,4),truncated(Normal(1.0,0.01),0.5,2)]
 
 bayesian_result = stan_inference(prob1,t,data,priors;num_samples=100,num_warmup=500,vars =(DiffEqBayes.StanODEData(),InverseGamma(4,1)))
-sdf  = CmdStan.read_summary(bayesian_result.model)
-@test sdf[sdf.parameters .== :theta1, :mean][1] ≈ 1.5 atol=1e-1
-@test sdf[sdf.parameters .== :theta2, :mean][1] ≈ 1.0 atol=1e-1
-@test sdf[sdf.parameters .== :theta3, :mean][1] ≈ 3.0 atol=1e-1
-@test sdf[sdf.parameters .== :theta4, :mean][1] ≈ 1.0 atol=1e-1
+@test mean(get(bayesian_result.chains,:theta_1)[1]) ≈ 1.5 atol=1e-1
+@test mean(get(bayesian_result.chains,:theta_2)[1]) ≈ 1.0 atol=1e-1
+@test mean(get(bayesian_result.chains,:theta_3)[1]) ≈ 3.0 atol=1e-1
+@test mean(get(bayesian_result.chains,:theta_4)[1]) ≈ 1.0 atol=1e-1

From 30b58f21cf680c1493ea588aa9d7198a28ff7a40 Mon Sep 17 00:00:00 2001
From: Vaibhav Dixit <vaibhavyashdixit@gmail.com>
Date: Sun, 11 Jul 2021 23:18:32 +0530
Subject: [PATCH 5/6] Update tests

---
 src/dynamichmc_inference.jl | 2 +-
 test/stan_string.jl         | 2 +-
 2 files changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/dynamichmc_inference.jl b/src/dynamichmc_inference.jl
index 834f1d4d..4aaf0623 100644
--- a/src/dynamichmc_inference.jl
+++ b/src/dynamichmc_inference.jl
@@ -103,5 +103,5 @@ function dynamichmc_inference(problem::DiffEqBase.DEProblem, algorithm, t, data,
     ℓ = TransformedLogDensity(trans, P)
     ∇ℓ = LogDensityProblems.ADgradient(AD_gradient_kind, ℓ)
     results = mcmc_with_warmup(rng, ∇ℓ, num_samples; mcmc_kwargs...)
-    merge((posterior = transform.(Ref(trans), results.chain), ), results)
+    merge((posterior = TransformVariables.transform.(Ref(trans), results.chain), ), results)
 end
diff --git a/test/stan_string.jl b/test/stan_string.jl
index e26d962e..e115dc49 100644
--- a/test/stan_string.jl
+++ b/test/stan_string.jl
@@ -1,4 +1,4 @@
-using CmdStan, Distributions, DiffEqBayes
+using Distributions, DiffEqBayes
 
 println("Starting the test")
 @test stan_string(Bernoulli(1)) == "bernoulli(1.0)"

From 9ed7f809ec4b60f5dd9a9eb55e85998c8ebf7cf9 Mon Sep 17 00:00:00 2001
From: Vaibhav Dixit <vaibhavyashdixit@gmail.com>
Date: Sun, 11 Jul 2021 23:31:00 +0530
Subject: [PATCH 6/6] Unexport stan_string

---
 test/stan_string.jl | 100 ++++++++++++++++++++++----------------------
 1 file changed, 50 insertions(+), 50 deletions(-)

diff --git a/test/stan_string.jl b/test/stan_string.jl
index e115dc49..0c52ed91 100644
--- a/test/stan_string.jl
+++ b/test/stan_string.jl
@@ -1,77 +1,77 @@
 using Distributions, DiffEqBayes
 
 println("Starting the test")
-@test stan_string(Bernoulli(1)) == "bernoulli(1.0)"
-@test stan_string(Bernoulli) == "bernoulli"
+@test DiffEqBayes.stan_string(Bernoulli(1)) == "bernoulli(1.0)"
+@test DiffEqBayes.stan_string(Bernoulli) == "bernoulli"
 
-@test stan_string(Binomial(5,0.3)) == "binomial(5, 0.3)"
-@test stan_string(Binomial) == "binomial"
+@test DiffEqBayes.stan_string(Binomial(5,0.3)) == "binomial(5, 0.3)"
+@test DiffEqBayes.stan_string(Binomial) == "binomial"
 
-@test stan_string(BetaBinomial(5,1,1)) == "beta_binomial(5, 1.0, 1.0)"
-@test stan_string(BetaBinomial) == "beta_binomial"
+@test DiffEqBayes.stan_string(BetaBinomial(5,1,1)) == "beta_binomial(5, 1.0, 1.0)"
+@test DiffEqBayes.stan_string(BetaBinomial) == "beta_binomial"
 
-@test stan_string(Hypergeometric(5,5,3)) == "hypergeometric(5, 5, 3)"
-@test stan_string(Hypergeometric) == "hypergeometric"
+@test DiffEqBayes.stan_string(Hypergeometric(5,5,3)) == "hypergeometric(5, 5, 3)"
+@test DiffEqBayes.stan_string(Hypergeometric) == "hypergeometric"
 
-@test stan_string(NegativeBinomial(5,0.3)) == "neg_binomial(5.0, 0.3)"
-@test stan_string(NegativeBinomial) == "neg_binomial"
+@test DiffEqBayes.stan_string(NegativeBinomial(5,0.3)) == "neg_binomial(5.0, 0.3)"
+@test DiffEqBayes.stan_string(NegativeBinomial) == "neg_binomial"
 
-@test stan_string(Poisson(5)) == "poisson(5.0)"
-@test stan_string(Poisson) == "poisson"
+@test DiffEqBayes.stan_string(Poisson(5)) == "poisson(5.0)"
+@test DiffEqBayes.stan_string(Poisson) == "poisson"
 
-@test stan_string(Normal(0,1)) == "normal(0.0, 1.0)"
-@test stan_string(Normal) == "normal"
+@test DiffEqBayes.stan_string(Normal(0,1)) == "normal(0.0, 1.0)"
+@test DiffEqBayes.stan_string(Normal) == "normal"
 
-@test stan_string(TDist(5)) == "student_t(5.0,0,1)"
-@test stan_string(TDist) == "student_t"
+@test DiffEqBayes.stan_string(TDist(5)) == "student_t(5.0,0,1)"
+@test DiffEqBayes.stan_string(TDist) == "student_t"
 
-@test stan_string(Cauchy(0,1)) == "cauchy(0.0, 1.0)"
-@test stan_string(Cauchy) == "cauchy"
+@test DiffEqBayes.stan_string(Cauchy(0,1)) == "cauchy(0.0, 1.0)"
+@test DiffEqBayes.stan_string(Cauchy) == "cauchy"
 
-@test stan_string(Laplace(0,1)) == "double_exponential(0.0, 1.0)"
-@test stan_string(Laplace) == "double_exponential"
+@test DiffEqBayes.stan_string(Laplace(0,1)) == "double_exponential(0.0, 1.0)"
+@test DiffEqBayes.stan_string(Laplace) == "double_exponential"
 
-@test stan_string(Distributions.Logistic(0,1)) == "logistic(0.0, 1.0)"
-@test stan_string(Distributions.Logistic) == "logistic"
+@test DiffEqBayes.stan_string(Distributions.Logistic(0,1)) == "logistic(0.0, 1.0)"
+@test DiffEqBayes.stan_string(Distributions.Logistic) == "logistic"
 
-@test stan_string(Gumbel(0,1)) == "gumbel(0.0, 1.0)"
-@test stan_string(Gumbel) == "gumbel"
+@test DiffEqBayes.stan_string(Gumbel(0,1)) == "gumbel(0.0, 1.0)"
+@test DiffEqBayes.stan_string(Gumbel) == "gumbel"
 
-@test stan_string(LogNormal(0,1)) == "lognormal(0.0, 1.0)"
-@test stan_string(LogNormal) == "lognormal"
+@test DiffEqBayes.stan_string(LogNormal(0,1)) == "lognormal(0.0, 1.0)"
+@test DiffEqBayes.stan_string(LogNormal) == "lognormal"
 
-@test stan_string(Chisq(5)) == "chi_square(5.0)"
-@test stan_string(Chisq) == "chi_square"
+@test DiffEqBayes.stan_string(Chisq(5)) == "chi_square(5.0)"
+@test DiffEqBayes.stan_string(Chisq) == "chi_square"
 
-@test stan_string(Exponential(5)) == "exponential(5.0)"
-@test stan_string(Exponential) == "exponential"
+@test DiffEqBayes.stan_string(Exponential(5)) == "exponential(5.0)"
+@test DiffEqBayes.stan_string(Exponential) == "exponential"
 
-@test stan_string(Gamma(2,3)) == "gamma(2.0, 3.0)"
-@test stan_string(Gamma) == "gamma"
+@test DiffEqBayes.stan_string(Gamma(2,3)) == "gamma(2.0, 3.0)"
+@test DiffEqBayes.stan_string(Gamma) == "gamma"
 
-@test stan_string(InverseGamma(2,3)) == "inv_gamma(2.0, 3.0)"
-@test stan_string(InverseGamma) == "inv_gamma"
+@test DiffEqBayes.stan_string(InverseGamma(2,3)) == "inv_gamma(2.0, 3.0)"
+@test DiffEqBayes.stan_string(InverseGamma) == "inv_gamma"
 
-@test stan_string(Weibull(1,1)) == "weibull(1.0, 1.0)"
-@test stan_string(Weibull) == "weibull"
+@test DiffEqBayes.stan_string(Weibull(1,1)) == "weibull(1.0, 1.0)"
+@test DiffEqBayes.stan_string(Weibull) == "weibull"
 
-@test stan_string(Frechet(1,1)) =="frechet(1.0, 1.0)"
-@test stan_string(Frechet) == "frechet"
+@test DiffEqBayes.stan_string(Frechet(1,1)) =="frechet(1.0, 1.0)"
+@test DiffEqBayes.stan_string(Frechet) == "frechet"
 
-@test stan_string(Rayleigh(5)) == "rayleigh(5.0)"
-@test stan_string(Rayleigh) == "rayleigh"
+@test DiffEqBayes.stan_string(Rayleigh(5)) == "rayleigh(5.0)"
+@test DiffEqBayes.stan_string(Rayleigh) == "rayleigh"
 
-@test stan_string(Pareto(2,3)) == "pareto(2.0, 3.0)"
-@test stan_string(Pareto) == "pareto"
+@test DiffEqBayes.stan_string(Pareto(2,3)) == "pareto(2.0, 3.0)"
+@test DiffEqBayes.stan_string(Pareto) == "pareto"
 
-@test stan_string(GeneralizedPareto(0,1,2)) == "pareto_type_2(0.0, 1.0, 2.0)"
-@test stan_string(GeneralizedPareto) == "pareto_type_2"
+@test DiffEqBayes.stan_string(GeneralizedPareto(0,1,2)) == "pareto_type_2(0.0, 1.0, 2.0)"
+@test DiffEqBayes.stan_string(GeneralizedPareto) == "pareto_type_2"
 
-@test stan_string(Beta(3,3)) == "beta(3.0, 3.0)"
-@test stan_string(Beta) == "beta"
+@test DiffEqBayes.stan_string(Beta(3,3)) == "beta(3.0, 3.0)"
+@test DiffEqBayes.stan_string(Beta) == "beta"
 
-@test stan_string(Uniform(1,4)) =="uniform(1.0, 4.0)"
-@test stan_string(Uniform) == "uniform"
+@test DiffEqBayes.stan_string(Uniform(1,4)) =="uniform(1.0, 4.0)"
+@test DiffEqBayes.stan_string(Uniform) == "uniform"
 
-@test stan_string(VonMises(0,2)) == "von_mises(0.0, 2.0)"
-@test stan_string(VonMises) == "von_mises"
+@test DiffEqBayes.stan_string(VonMises(0,2)) == "von_mises(0.0, 2.0)"
+@test DiffEqBayes.stan_string(VonMises) == "von_mises"