使用 deSolve package- ddply 函数模拟数据

Question

我无法在 deSolve 包中为每个受试者应用微分方程求解器 ID 来计算 2 室静脉输液模型中的药物量。我能够设置代码（如下），因此它解决了一个主题。我需要有关如何使用 ddply 包将它应用于数据框中的每个主题的帮助。

下面是一个名为 simeventdfi 的数据框，其中包含计算所需的剂量事件。

library(deSolve)

第一步：为每个主题提供事件 df

 simeventdfi <- 
 ID  var time value method
 1    3   0.0   6   rep
 1    3  16.7   0   rep
 1    4  0.0    2.4 rep
 1    4  16.7   2.4 rep
 1    5  0.0    10  rep
 1    5  16.7   10  rep
 1    6  0.0    1   rep
 1    6  16.7   1   rep
 1    7  0.0    25  rep
 1    7  16.7   25  rep
 2    3  0.0    6   rep
 2    3  16.7   0   rep
 2    4  0.0    2.4 rep
 2    4  16.7   2.4 rep
 2    5  0.0    10  rep
 2    5  16.7   10  rep
 2    6  0.0    1   rep
 2    6  16.7   1   rep
 2    7  0.0    25  rep
 2    7  16.7   25  rep

第 2 步：指定模拟值的时间

   simtimes <- sort(unique(c(simeventdfi$time, seq(0,10,1))))

第 3 步：包含隔间 (A) 中量的微分方程的函数

  #THETAin is empty for this model   
   simthetai <- NULL

 DES <- function(T, A, THETAin)
  {

    RATE  <- A[3]    #Dose rate (time dependent)
    CL <- A[4]         #Time dependent                           
    V1 <- A[5]         #Time dependent
    Q  <- A[6]         #Time dependent
    V2 <- A[7]         #Time dependent

    dA1 <- RATE -Q/V1*A[1] +Q/V2*A[2] -CL/V1*A[1]  #Central compartment 
    dA2 <- Q/V1*A[1] - Q/V2*A[2]     #Peripheral compartment 

    RATE <- 0  #Set rate to zero so doesn't change unless event
    dCL <- 0 
    dV1 <- 0
    dQ <- 0
    dV2 <- 0   
   list(c(dA1,dA2,RATE,dCL,dV1,dQ,dV2))   #List of derivatives
     }

第 4 步：由于我能够将函数代码设置为仅解决一个主题 ID，因此我必须从上面给出的 simeventdfi 中提取一个主题 ID。但是，我需要有关如何使用 ddply 为每个主题 ID 应用微分方程求解器的帮助。现在，我将为一个主题对上面的 simeventdfi 进行子集化，以演示该函数的工作原理。

simeventdfi <- subset(simeventdfi,ID==1)

设置初始值 - 隔间和时间相关参数

    A_0i <- c("A1"=0,"A2"=0,
      "Rate"=simeventdfi$value[1],
      "CL"=simeventdfi$value[simeventdfi$var==4 & simeventdfi$time==0],
      "V1"=simeventdfi$value[simeventdfi$var==5 & simeventdfi$time==0],
      "Q"= simeventdfi$value[simeventdfi$var==6 & simeventdfi$time==0],
      "V2"=simeventdfi$value[simeventdfi$var==7 & simeventdfi$time==0]) 
   print(A_0i)

第 5 步：运行微分方程求解器并在数据框中获取结果

    simdatadfi <- as.data.frame(ode(A_0i, simtimes, DES, simthetai, events=list(data=simeventdfi), method="lsoda"))

它只适用于一个主题。我需要有关如何为 simeventdfi.

中的每个主题应用步骤 5 中的微分方程求解器的帮助

Answer 1

您可以使用 dplyr。只需使用 group_by(ID) 为每个主题应用微分方程。显然 ode 不喜欢 dplyr 数据格式，因此您需要在 do 调用中使用 data.frame(.) 而不仅仅是 .。

require(dplyr)
simeventdfi %>% 
  group_by(ID) %>% 
  do(
    ode(A_0i, 
        simtimes, 
        DES, 
        simthetai, 
        events=list(data=data.frame(.)), 
        method="lsoda") %>%
      data.frame
  ) %>% 
  data.frame

Answer 2

要将微分方程求解器应用于每个主题：

首先：将第4步写成一个函数：

  simulate.conc <- function(simeventdfi) { 

      #Initial values - compartments and time-dependent parameters
      A_0i <- c("A1"=0,"A2"=0,
        "Rate"=simeventdfi$value[1],
        "CL"=simeventdfi$value[simeventdfi$var==4 & simeventdfi$time==0],
        "V1"=simeventdfi$value[simeventdfi$var==5 & simeventdfi$time==0],
        "Q"= simeventdfi$value[simeventdfi$var==6 & simeventdfi$time==0],
        "V2"=simeventdfi$value[simeventdfi$var==7 & simeventdfi$time==0]) 

   #Run differential equation solver 
     simdata <- ode(A_0i, simtimes, DES, simthetai, events=list(data=simeventdfi), method="lsoda")

 #Process the simulated output 
  simdata <- as.data.frame(simdata)

  #Concentration in the central compartment
  simdata$conc <- simdata$A1/simdata$V1    #Amount to concentration

  #remove any duplicated datapoints
  simdata <- simdata[!duplicated(simdata), ]
 }

其次：使用 ddply

将函数应用于 simeventdfi 中的每个 ID

 library(plyr)
 simdata_DES <- ddply(simeventdfi, .(ID), simulate.conc)

简单！

使用 deSolve package- ddply 函数模拟数据

Simulating data using deSolve package- ddply function

r

plyr

differential-equations