added conditional dependence of t on s. Converging but not well. Need to setup without it and compare estimates.

check_reversion
Will King 3 years ago
parent 611fd48841
commit 8a8c49f9ca

@ -23,9 +23,6 @@ A,B = f(μ,σ)
(s,t) ~ (bernoulli(p), Gamma(A,B)) (s,t) ~ (bernoulli(p), Gamma(A,B))
""" """
# ╔═╡ d3c53cb5-dd7a-4da0-82d6-52611a1a9e1c
# ╔═╡ fb1ce990-915e-46eb-88a6-ac3c2bc55f4e # ╔═╡ fb1ce990-915e-46eb-88a6-ac3c2bc55f4e
begin #test params begin #test params
n = 100 n = 100
@ -54,25 +51,25 @@ begin
X = hcat(τ, x₁, x₂, x₃) X = hcat(τ, x₁, x₂, x₃)
end end
# ╔═╡ 8735884f-efdc-4421-a205-06b68946b564
begin #build statuses
ξ = [-5 6 -7 3 -2 4]'
Ξ = logistic.(X * ξ)
s = rand.(Bernoulli.(Ξ))
end
# ╔═╡ 7cee4dad-17b4-4cf8-8fa7-631722c4f9bf # ╔═╡ 7cee4dad-17b4-4cf8-8fa7-631722c4f9bf
begin #Build Durations begin #Build Durations
β = [1 2 3 4 5 1]' β = [1 2 3 4 5 1]'
ν = [6 5 4 3 2 1]' ν = [6 5 4 3 2 1]'
α = [-1,1]
μ = X * β
σ = X * ν μ = X * β .+ [α[x+1] for x in s]
η = (X * ν).^2
scale = σ.^2 ./ μ
shape = μ ./ scale
t = rand.(Gamma.(shape,scale))
end
# ╔═╡ 8735884f-efdc-4421-a205-06b68946b564 t = rand.(Normal.(μ,η))
begin #build statuses
ξ = [-5 6 -7 3 -2 4]'
Ξ = logistic.(X * ξ)
s = rand.(Bernoulli.(Ξ))
end end
# ╔═╡ f0687986-ac5f-47c7-ad79-07318d770d17 # ╔═╡ f0687986-ac5f-47c7-ad79-07318d770d17
@ -83,7 +80,7 @@ plot(Ξ,s,seriestype=:scatter, legend=nothing)
n,k=size(X) n,k=size(X)
#model S #model S
ξ ~ MvNormal(zeros(k),1) ξ ~ MvNormal(zeros(k),2)
σₛ ~ Exponential(1) σₛ ~ Exponential(1)
μ = X*ξ μ = X*ξ
@ -92,14 +89,25 @@ plot(Ξ,s,seriestype=:scatter, legend=nothing)
s .~ Bernoulli.(p) s .~ Bernoulli.(p)
#model t #model t
β ~ MvNormal(zeros(k),2)
ν ~ MvNormal(zeros(k),2)
α ~ MvNormal([0,0], 2)
μ = X * β .+ [α[x+1] for x in s]
η = (X * ν).^2
t .~ Normal.(μ,η)
end end
# ╔═╡ a9362616-b992-4100-b000-48c8ae7b4d4b
a =[1,2]
# ╔═╡ edad9362-fbd3-445a-857d-9687c1adbddd # ╔═╡ edad9362-fbd3-445a-857d-9687c1adbddd
testmod = MultiVariate(s,t,X) testmod = MultiVariate(s,t,X)
# ╔═╡ d30d54df-260e-4826-89e3-de5434f8a730 # ╔═╡ d30d54df-260e-4826-89e3-de5434f8a730
chain = sample(testmod,NUTS(),200) chain = sample(testmod,Gibbs(NUTS(),PG(20,:α)),1000)
# ╔═╡ 846a1840-5734-4ed3-9d9c-01f063cd0057 # ╔═╡ 846a1840-5734-4ed3-9d9c-01f063cd0057
plot(chain) plot(chain)
@ -1732,9 +1740,8 @@ version = "1.4.1+0"
""" """
# ╔═╡ Cell order: # ╔═╡ Cell order:
# ╠═c39be064-648a-11ed-0677-ebe624104c69 # ╟─c39be064-648a-11ed-0677-ebe624104c69
# ╠═add4bb20-1da8-45e7-a40f-1b41c149d5f0 # ╠═add4bb20-1da8-45e7-a40f-1b41c149d5f0
# ╠═d3c53cb5-dd7a-4da0-82d6-52611a1a9e1c
# ╠═fb1ce990-915e-46eb-88a6-ac3c2bc55f4e # ╠═fb1ce990-915e-46eb-88a6-ac3c2bc55f4e
# ╠═d75ac31e-9524-4ad7-b350-0a6dabdf22d2 # ╠═d75ac31e-9524-4ad7-b350-0a6dabdf22d2
# ╟─a84568ec-f443-4b36-83ba-baf5aceae209 # ╟─a84568ec-f443-4b36-83ba-baf5aceae209
@ -1742,6 +1749,7 @@ version = "1.4.1+0"
# ╠═8735884f-efdc-4421-a205-06b68946b564 # ╠═8735884f-efdc-4421-a205-06b68946b564
# ╠═f0687986-ac5f-47c7-ad79-07318d770d17 # ╠═f0687986-ac5f-47c7-ad79-07318d770d17
# ╠═fd7264ce-00af-49a3-91c9-2277bad2b3e9 # ╠═fd7264ce-00af-49a3-91c9-2277bad2b3e9
# ╠═a9362616-b992-4100-b000-48c8ae7b4d4b
# ╠═edad9362-fbd3-445a-857d-9687c1adbddd # ╠═edad9362-fbd3-445a-857d-9687c1adbddd
# ╠═d30d54df-260e-4826-89e3-de5434f8a730 # ╠═d30d54df-260e-4826-89e3-de5434f8a730
# ╠═846a1840-5734-4ed3-9d9c-01f063cd0057 # ╠═846a1840-5734-4ed3-9d9c-01f063cd0057

Loading…
Cancel
Save