diff --git a/Julia-turing models/Model.jl b/Julia-turing models/Model.jl
new file mode 100644
index 0000000..60d94d7
--- /dev/null
+++ b/Julia-turing models/Model.jl	
@@ -0,0 +1,62 @@
+
+using Turing, StatsFuns
+using Distributions, StatsPlots
+
+using FillArrays
+
+begin #creating synthetic data
+	βs = [1 2 3]'
+	x = Matrix([1:10 1:2:20 1:3:30])
+	X = (x .- mean(x,dims=1)) ./ std(x,dims=1)
+	t = [x for x=1:0.5:5.5]
+	s = rand([1,2],10)
+	σ= [2 3]
+	rand_draw1 = randn(10)
+
+	y = x*βs .+ σ[s]
+end
+
+@model function JointDurationStateModel(
+	DeviationFromExpectedDuration,
+	ConclusionStatus,
+	SnapshotState,
+	CurrentDuration,
+)
+	# get dimensions
+	n,k = size(SnapshotState)
+	
+	#hyperpriors priors
+
+	#Heirarchal parameters
+	#β ~ MvNormal(Fill(0,k),2)
+	η ~ MvNormal(Fill(0,k),2)
+
+	#Direct Priors
+	#σ_DFED ~ filldist(Exponential(1),2) #TODO: check implication of this form
+	
+	#model
+	#μ = SnapshotState * β
+	p = StatsFuns.logistic.(SnapshotState * η)
+	
+	#estimate ConclusionStatus model
+	ConclusionStatus .~ Bernoulli(p)
+	#Estimate DFED model
+	#=
+	for i in eachindex(ConclusionStatus)
+		DeviationFromExpectedDuration ~ Normal(
+			μ[ConclusionStatus[i]],
+			σ_DFED[ConclusionStatus[i]]
+		)
+	end
+	=#
+end
+
+model = JointDurationStateModel(y,s,X,t)
+prior = JointDurationStateModel(fill(missing,size(y)),fill(missing,size(s)),X,t)
+
+chain = sample(model,NUTS(0.85),2000)
+
+
+
+plot(chain)
+