首頁|安裝|文檔|教程|開發商
rodeo是一個快速的 Python 函式庫,它使用機率數值來求解常微分方程 (ODE)。也就是說,大多數 ODE 解算器(例如歐拉法)會在步長網格上產生 ODE 的確定性近似值
rodeo為許多機率求解器常見的非線性貝葉斯濾波範例提供了一個輕量級且可擴展的近似系列(Tronarp 等人(2018))。首先將高斯過程置於 ODE 解的前面,並在求解器逐步遍歷網格時按順序更新它。 rodeo基於jax構建,允許即時編譯和自動微分。 jax的 API 幾乎等同於numpy的 API。
rodeo提供了兩種主要工具:一種用於近似 ODE 解,另一種用於參數推斷。對於前者我們提供:
solve
:使用非線性貝葉斯濾波範例的機率 ODE 求解器的實現。對於後者,我們提供似然近似方法:
basic
:基本似然近似方法的實現(詳細資訊可在 Wu 和 Lysy (2023) 中找到)。fenrir
:Fenrir 的實作(Tronarp 等人(2022))。marginal_mcmc
:Chkrebtii 方法的 MCMC 實作(Chkrebtii 等人(2016))。dalton
:資料自適應 ODE 似然近似的實現(Wu 和 Lysy (2023))。其使用的詳細範例可以在文件部分找到。請注意,這是rodeo的僅jax版本。使用各種其他後端的舊版本,請參閱此處。
從 GitHub 下載儲存庫,然後使用setup.cfg
腳本進行安裝:
git clone https://github.com/mlysy/rodeo.git
cd rodeo
pip install .
請先前往 readthedocs 查看以下範例的渲染文件。
有關解決簡單 ODE 問題的快速入門教學。
求解高階 ODE 的範例。
求解困難的混沌 ODE 的範例。
我們使用拉普拉斯近似的參數推斷問題的範例。
在本演練中,我們展示如何使用機率求解器來求解 ODE 並進行參數推論。我們將首先說明求解 ODE 的設定。為此,我們考慮以下一階 ODE 範例( FitzHugh-Nagumo模型),
解決方案在哪裡
根據 (Wu and Lysy (2023)) 的記號,我們有
為了
import jax
import jax . numpy as jnp
import rodeo
def fitz_fun ( X , t , ** params ):
"FitzHugh-Nagumo ODE in rodeo format."
a , b , c = params [ "theta" ]
V , R = X [:, 0 ]
return jnp . array (
[[ c * ( V - V * V * V / 3 + R )],
[ - 1 / c * ( V - a + b * R )]]
)
def fitz_init ( x0 , theta ):
"FitzHugh-Nagumo initial values in rodeo format."
x0 = x0 [:, None ]
return jnp . hstack ([
x0 ,
fitz_fun ( X = x0 , t = 0. , theta = theta ),
jnp . zeros_like ( x0 )
])
W = jnp . array ([[[ 0. , 1. , 0. ]], [[ 0. , 1. , 0. ]]]) # LHS matrix of ODE
x0 = jnp . array ([ - 1. , 1. ]) # initial value for the ODE-IVP
theta = jnp . array ([ .2 , .2 , 3 ]) # ODE parameters
X0 = fitz_init ( x0 , theta ) # initial value in rodeo format
# Time interval on which a solution is sought.
t_min = 0.
t_max = 40.
# --- Define the prior process -------------------------------------------
n_vars = 2 # number of variables in the ODE
n_deriv = 3 # max number of derivatives
sigma = jnp . array ([ .1 ] * n_vars ) # IBM process scale factor
# --- data simulation ------------------------------------------------------
n_steps = 800 # number of evaluations steps
dt = ( t_max - t_min ) / n_steps # step size
# generate the Kalman parameters corresponding to the prior
prior_Q , prior_R = rodeo . prior . ibm_init (
dt = dt_sim ,
n_deriv = n_deriv ,
sigma = sigma
)
# Produce a Pseudo-RNG key
key = jax . random . PRNGKey ( 0 )
Xt , _ = rodeo . solve_mv (
key = key ,
# define ode
ode_fun = fitz_fun ,
ode_weight = W ,
ode_init = X0 ,
t_min = t_min ,
t_max = t_max ,
theta = theta , # ODE parameters added here
# solver parameters
n_steps = n_steps ,
interrogate = rodeo . interrogate . interrogate_kramer ,
prior_weight = prior_Q ,
prior_var = prior_R
)
我們將求解器的解與scipy庫中odeint
提供的確定性解進行比較。
我們也提供了一個求解高階 ODE 和混沌 ODE 的範例。
我們現在轉向參數推斷問題。 rodeo包含描述部分中總結的幾種似然近似方法。在這裡,我們將使用basic
似然近似方法。假設透過模型模擬觀察結果
在哪裡basic
似然近似
def fitz_logprior ( upars ):
"Logprior on unconstrained model parameters."
n_theta = 5 # number of ODE + IV parameters
lpi = jax . scipy . stats . norm . logpdf (
x = upars [: n_theta ],
loc = 0. ,
scale = 10.
)
return jnp . sum ( lpi )
def fitz_loglik ( obs_data , ode_data , ** params ):
"""
Loglikelihood for measurement model.
Args:
obs_data (ndarray(n_obs, n_vars)): Observations data.
ode_data (ndarray(n_obs, n_vars, n_deriv)): ODE solution.
"""
ll = jax . scipy . stats . norm . logpdf (
x = obs_data ,
loc = ode_data [:, :, 0 ],
scale = 0.005
)
return jnp . sum ( ll )
def constrain_pars ( upars , dt ):
"""
Convert unconstrained optimization parameters into rodeo inputs.
Args:
upars : Parameters vector on unconstrainted scale.
dt : Discretization grid size.
Returns:
tuple with elements:
- theta : ODE parameters.
- X0 : Initial values in rodeo format.
- Q, R : Prior matrices.
"""
theta = jnp . exp ( upars [: 3 ])
x0 = upars [ 3 : 5 ]
X0 = fitz_init ( x0 , theta )
sigma = upars [ 5 :]
Q , R = rodeo . prior . ibm_init (
dt = dt ,
n_deriv = n_deriv ,
sigma = sigma
)
return theta , X0 , Q , R
def neglogpost_basic ( upars ):
"Negative logposterior for basic approximation."
# solve ODE
theta , X0 , prior_Q , prior_R = constrain_pars ( upars , dt_sim )
# basic loglikelihood
ll = rodeo . inference . basic (
key = key ,
# ode specification
ode_fun = fitz_fun ,
ode_weight = W ,
ode_init = X0 ,
t_min = t_min ,
t_max = t_max ,
theta = theta ,
# solver parameters
n_steps = n_steps ,
interrogate = rodeo . interrogate . interrogate_kramer ,
prior_weight = prior_Q ,
prior_var = prior_R ,
# observations
obs_data = obs_data ,
obs_times = obs_times ,
obs_loglik = fitz_loglik
)
return - ( ll + fitz_logprior ( upars ))
這是示範用法的基本範例。我們建議使用更複雜的似然近似,將解的不確定性傳播到似然近似,例如fenrir
、 marginal_mcmc
和dalton
。更多詳情請參考參數推斷教學。
以下是/examples/
中在rodeo中發現的各種似然近似產生的一些結果:
可以透過以下命令運行單元測試:
cd tests
python -m unittest discover -v
或者,安裝tox ,然後從rodeo
輸入命令列: tox
。
HTML 文件可以從根資料夾編譯:
pip install .[docs]
cd docs
make html
這將在docs/build
中建立文件。