0

生物学的プロセスに関連する一連の方程式を解こうとしています。1 つの式 (約 5 つ) は、 の形の薬物動態 (PK) 曲線用ですC = Co(exp(k1*t)-exp(k2*t)。この方程式の導関数を、いくつかの酵素結合方程式と期待どおりではない初期結果とともに同時に解く必要があります。トラブルシューティングの後、desolve ode 関数を使用して k が負の場合、PK 導関数はそれ自体では数値的に積分されないことに気付きました。ode 関数ですべてのメソッド (lsode、lsoda など) を試しましたが、成功しませんでした。rtol を調整してみましたが、解決しません。

私が調査すべき deSolve ode 関数に代わるものはありますか? または、この問題を解決する別の方法はありますか?

以下は、問題を説明するための単純化された方程式を含むコードです。k が負の場合、積分解は解析結果と一致しません。k が正の場合、結果は期待どおりです。

最初の画像、k=0.2 の結果: k が正の場合、解析結果と統合結果は一致します k が正の場合、解析結果と統合結果は一致します

2 番目の画像、k=-0.2 の結果: k が負の場合、統合された結果は分析と一致しません k が負の場合、積分結果は解析結果と一致しません

library(deSolve)

abi <- function(t, state, parameters) {
    with(as.list(c(state, parameters)), {
        
        dI <- k*exp(k*t)
                list(c(dI))
    })
}

k <- c(-0.2)

times <- seq(0, 24, by = 1)

I_analytical <- exp(k*times)

parameters <- c(k)
state <- c(I = 0)

out <- ode(y = state, times = times, func = abi, parms = parameters)

plot(out)
points(I_analytical ~ times)

初期条件は上記の例を簡単に解決することが指摘されており、非常に役に立ちます。これは、正確に積分できない方程式です。いくつかの異なる初期条件を試しましたが、実際には成功しませんでした。

library(deSolve)

## Chaos in the atmosphere
CYP <- function(t, state, parameters) {
    with(as.list(c(state, parameters)), {
        #dE <- ksyn - (kdeg * E) + (k2 * EI) - (k1 * E * I)
        #dEI <- (k1 * E * I) - (k2 * EI) + (k4 * EIstar) - (k3 * EI)
        #dEIstar <- (k3 * EI) - (k4 * EIstar)
        #dOcc <- dEI + dEIstar
        dI <- a*tau1*exp(tau1*t) + b*tau2*exp(tau2*t) + c*tau3*exp(tau3*t)
            #list(c(dE, dEI, dEIstar, dOcc, dI))
          list(c(dI))
    })
}

ifit <- c(-0.956144311,0.82619445,0.024520276,-0.913499862,-0.407478829,-0.037174745)
a = ifit[1]
b = ifit[2]
c = ifit[3]
tau1 = ifit[4]
tau2 = ifit[5]
tau3 = ifit[6]


parameters <- c(ksyn = 0.82, kdeg = 0.02, k1 = 2808, k2 = 370.66, k3 = 2.12, k4 = 0.017, a, b, c, tau1, tau2, tau3)

#state <- c(E = 41, EI = 0, EIstar = 0, Occupancy = 0, I = 0.0)
state <- c(I=-0.01)
times <- seq(0, 24, by = .1)
out <- ode(y = state, times = times, func = CYP, parms = parameters)

I_analytical <- a*exp(tau1*times) + b*exp(tau2*times) + c*exp(tau3*times)

plot(out)
points(I_analytical ~ times)

ターゲット カーブとオード ソリューション ライン。

4

2 に答える 2

0

初期値は

state <- c(I= a + b + c)
#state <- c(I = 1)
于 2021-04-30T15:23:52.213 に答える