如何在 Julia 中调整贝叶斯 ODE 的超参数?
How to tune the hyperparameters of a Bayesian ODE fit in Julia?
我一直在尝试使用不同的 ODE 方程来复制 https://diffeqflux.sciml.ai/dev/examples/BayesianNODE_NUTS/,但我收到的结果没有不确定性量化,是不是因为我的初始值 u0 更高:
你能告诉我哪里出了问题吗?
using DiffEqFlux, OrdinaryDiffEq, Flux, Optim, Plots, AdvancedHMC, MCMCChains
using JLD, StatsPlots
function Arps!(du,u,p,t)
y= u[1]
#x, y = u
# Di,b,n,tau = p
n,tau = p
#du[1]=dx=-(x * Di * x^b)
du[1]=dy=-(n *((t^n)/tau) * y/t)
end
tspan=(1.0,50.0)
tsteps = 1:1:50
u0 = [16382.9]
p=[0.48,15.92]
prob_trueode = ODEProblem(Arps!,u0,tspan,p)
ode_data = Array(solve(prob_trueode, Tsit5(), saveat = tsteps))
ode_data =ode_data[1,:]
dudt= FastChain(FastDense(1, 30, tanh),
FastDense(30, 1))
prob_neuralode = NeuralODE(dudt, tspan, Tsit5(), saveat = tsteps)
function predict_neuralode(p)
Array(prob_neuralode(u0, p))
end
function loss_neuralode(p)
pred = predict_neuralode(p)
loss = sum(abs2, ode_data .- pred)
return loss, pred
end
l(θ) = -sum(abs2, ode_data .- predict_neuralode(θ)) - sum(θ .* θ)
function dldθ(θ)
x,lambda = Flux.Zygote.pullback(l,θ)
grad = first(lambda(1))
return x, grad
end
metric = DiagEuclideanMetric(length(prob_neuralode.p))
h = Hamiltonian(metric, l, dldθ)
integrator = Leapfrog(find_good_stepsize(h, Float64.(prob_neuralode.p)))
prop = AdvancedHMC.NUTS{MultinomialTS, GeneralisedNoUTurn}(integrator)
adaptor = StanHMCAdaptor(MassMatrixAdaptor(metric), StepSizeAdaptor(0.45, prop.integrator))
samples, stats = sample(h, prop, Float64.(prob_neuralode.p), 500, adaptor, 500; progress=true)
losses = map(x-> x[1],[loss_neuralode(samples[i]) for i in 1:length(samples)])
################### RETRODICTED PLOTS: TIME SERIES #################
pl = scatter(tsteps, ode_data, color = :red, label = "Data: Var1", xlabel = "t", title = "Spiral Neural ODE")
for k in 1:300
resol = predict_neuralode(samples[100:end][rand(1:400)])
plot!(tsteps,resol[1,:], alpha=0.04, color = :red, label = "")
end
idx = findmin(losses)[2]
prediction = predict_neuralode(samples[idx])
plot!(tsteps,prediction[1,:], color = :black, w = 2, label = "")
最可能的原因是后验样本的损失函数幅度过高,因此后验样本结果超出范围并且在您的图上不可见。
这可以通过 (a) 添加神经 ODE 输出的比例因子并确保损失函数不是从非常高的幅度开始或 (b) 增加神经网络中的层数来解决架构/改变激活函数。
通过在Neural ODE中添加比例因子,我得到了如下图所示的良好结果:
我一直在尝试使用不同的 ODE 方程来复制 https://diffeqflux.sciml.ai/dev/examples/BayesianNODE_NUTS/,但我收到的结果没有不确定性量化,是不是因为我的初始值 u0 更高:
你能告诉我哪里出了问题吗?
using DiffEqFlux, OrdinaryDiffEq, Flux, Optim, Plots, AdvancedHMC, MCMCChains
using JLD, StatsPlots
function Arps!(du,u,p,t)
y= u[1]
#x, y = u
# Di,b,n,tau = p
n,tau = p
#du[1]=dx=-(x * Di * x^b)
du[1]=dy=-(n *((t^n)/tau) * y/t)
end
tspan=(1.0,50.0)
tsteps = 1:1:50
u0 = [16382.9]
p=[0.48,15.92]
prob_trueode = ODEProblem(Arps!,u0,tspan,p)
ode_data = Array(solve(prob_trueode, Tsit5(), saveat = tsteps))
ode_data =ode_data[1,:]
dudt= FastChain(FastDense(1, 30, tanh),
FastDense(30, 1))
prob_neuralode = NeuralODE(dudt, tspan, Tsit5(), saveat = tsteps)
function predict_neuralode(p)
Array(prob_neuralode(u0, p))
end
function loss_neuralode(p)
pred = predict_neuralode(p)
loss = sum(abs2, ode_data .- pred)
return loss, pred
end
l(θ) = -sum(abs2, ode_data .- predict_neuralode(θ)) - sum(θ .* θ)
function dldθ(θ)
x,lambda = Flux.Zygote.pullback(l,θ)
grad = first(lambda(1))
return x, grad
end
metric = DiagEuclideanMetric(length(prob_neuralode.p))
h = Hamiltonian(metric, l, dldθ)
integrator = Leapfrog(find_good_stepsize(h, Float64.(prob_neuralode.p)))
prop = AdvancedHMC.NUTS{MultinomialTS, GeneralisedNoUTurn}(integrator)
adaptor = StanHMCAdaptor(MassMatrixAdaptor(metric), StepSizeAdaptor(0.45, prop.integrator))
samples, stats = sample(h, prop, Float64.(prob_neuralode.p), 500, adaptor, 500; progress=true)
losses = map(x-> x[1],[loss_neuralode(samples[i]) for i in 1:length(samples)])
################### RETRODICTED PLOTS: TIME SERIES #################
pl = scatter(tsteps, ode_data, color = :red, label = "Data: Var1", xlabel = "t", title = "Spiral Neural ODE")
for k in 1:300
resol = predict_neuralode(samples[100:end][rand(1:400)])
plot!(tsteps,resol[1,:], alpha=0.04, color = :red, label = "")
end
idx = findmin(losses)[2]
prediction = predict_neuralode(samples[idx])
plot!(tsteps,prediction[1,:], color = :black, w = 2, label = "")
最可能的原因是后验样本的损失函数幅度过高,因此后验样本结果超出范围并且在您的图上不可见。
这可以通过 (a) 添加神经 ODE 输出的比例因子并确保损失函数不是从非常高的幅度开始或 (b) 增加神经网络中的层数来解决架构/改变激活函数。
通过在Neural ODE中添加比例因子,我得到了如下图所示的良好结果: