我有一个关于生物系统的假设，我需要使用程序“deSolve”来表达微分方程和参数，以便模拟该系统。但我不知道我的代码出了什么问题。
这是我正在使用的系统：

这是模拟图应该看起来像：

这是我的图表目前的样子：

这是我目前拥有的代码：

#1. Define Initial conditions
#2. Define time-steps
#3. Define differential-equations for all states
#4.Simulate the model

#----------------------------------------------------------------------------
library(deSolve)

#Define initial conditions
states <- c(R=1, 
            Rp=0,
            RS=1,
            RSp=0,
            S=1)

#Define time steps
times <- c(seq(0, 20, 0.1))

#Define ordinary differential equations for the model
model1 <- function(time, states, parameters) {
  with (as.list(c(states, parameters)), {
    dR = -(R*k1*S) + (Rp*k2) + (Rp*kfeed*RSp)      
    dRp = (R*k1*S) - (Rp*k2) - (Rp*kfeed*RSp)
    dRS = -(RS*k5*Rp) + (RSp*k4)
    dRSp = (RS*k5*Rp) - (RSp*k4) - (Rp*kfeed*RSp)
    dS = 0;
    return(list(c(dR, dRp, dRS, dRSp, dS)))
  })
}

#Simulate the model
simModel1 <- function(parameters) {
  return(as.data.frame(ode(y = states, times = times, func = model1, parms = parameters)))
}

pstart1 <- c(k1=1, 
             k2=0.001,
             kfeed=100,
             k4=0.01,
             k5=0.01)

sim1 <- simModel1(pstart1)

plot(sim1$time, sim1$Rp, type = "l", xlab = "time", ylab = "Rp", col = 'violet')