return mcmcchains

Project.toml

@@ -11,6 +11,7 @@ JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
 JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 MAT = "23992714-dd62-5051-b70f-ba57cb901cac"
+MCMCChains = "c7f686f2-ff18-58e9-bc7b-31028e88f75d"
 MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
 OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
 Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"

src/COBREXA.jl

@@ -13,6 +13,7 @@ using Random
 using Serialization
 using SparseArrays
 using Statistics
+using MCMCChains
 
 import Base: findfirst, getindex, show
 import Pkg

src/analysis/sampling/hit_and_run.jl

@@ -73,7 +73,7 @@ function hit_and_run(
     map(fetch, save_at.(workers, :ubs, Ref(:($ubs))))
 
     # do in parallel! 
-    chains = dpmap(
+    samples = dpmap(
         x -> :($COBREXA._serial_hit_and_run(ws, lbs, ubs, $samplesize, $keepevery, $N)),
         CachingPool(workers),
         1:nchains,
@@ -84,6 +84,13 @@ function hit_and_run(
     map(fetch, remove_from.(workers, :lbs))
     map(fetch, remove_from.(workers, :ubs))
 
+    # not sure how to do this better
+    vals = zeros(samplesize, length(lbs), nchains)
+    for c in 1:nchains
+        vals[:, :, c] = samples[c]'
+    end
+    chains = Chains(vals, reactions(model))
+
     return chains
 end