0

元の質問を間違って行ったので、こちらがより良いバージョンです。

R を使用してモデルを生成したいと思います。モデルの要点 --> ポリマーはさまざまな速度で成長または収縮できます。時折、成長率が一定に保たれているにもかかわらず、縮小率が 20 倍に増加します。これは「壊滅的な状態」と見なされます。状態は、特定の割合で「壊滅的な状態」との間で切り替わります。問題は、ポリマーの長さが時間に対してどのように変化するかということです。これが私が考えていることです：

初期化:

長さ 0 の 5 つのポリマー (列インデックスで示されます)

rnumbers <- data.frame(replicate(5,runif(20000, 0, 1)))

rstate <- rnumbers  # copy the structure
rstate[] <- NA      # preserve structure with NA's
# Init:
rstate[1, ] <- c(0)

> head(rstate)
  X1 X2 X3 X4 X5
1  0  0  0  0  0
2 NA NA NA NA NA
3 NA NA NA NA NA
4 NA NA NA NA NA
5 NA NA NA NA NA
6 NA NA NA NA NA

シミュレーションを 200 秒間実行したい

レートの設定:

dt <- c(.01) # how time is being divided
A <- dt*1 # probability growth will occur under normal conditions 
B <- dt*.25 # probability shrinking will occur under normal conditions
C <- dt*.03 # probability "catastrophic state" will occur
D <- dt*.3 # probability normal state will occur once in "catastrophic state"
E <- dt*5 # probability shrinking will occur "catastrophic state"

通常の状態では成長の確率が収縮の確率を上回っていますが、「壊滅的な状態」では収縮が優勢です。また、データ フレームの 20000 行の dt = .01 は、最大 200 になります

壊滅的な状態への切り替えを考慮しないと、コードは次のようになります。

library("Rcell")
rnumbers <- data.frame(replicate(5,runif(20000, 0, 1)))
dt <- c(.01)
A <- dt*1
B <- dt*.25
C <- dt*.03
D <- dt*.3
E <- dt*5

rstate <- rnumbers  # copy the structure
rstate[] <- NA      # preserve structure with NA's
# Init:
rstate[1, ] <- c(0)

step_generator <- function(col, rnum){
    for (i in 2:length(col) ){
            if( rnum[i] < B) { 
                col[i] <- -1
                }
                       else { 
                        if (rnum[i] < A) {
                            col[i] <- 1 
                            }
                              else {
                                col[i] <- 0 
                                }
                                }
                        }
    return(col)
    }
#  Run for each column index:
for(cl in 1:5){ rstate[ , cl] <- 
                        step_generator(rstate[,cl], rnumbers[,cl]) }

rstate1 <- transform(rstate, time = rep(dt))
rstate2 <- transform(rstate1, cumtime = cumsum(time))

cum_sum <- apply(rstate2, 2, cumsum)
cum_sum <- as.data.frame(cum_sum)

dev.new(width=5, height=5)  
cplot(cum_sum, X2 ~ time, geom = "line")

このコードを実行すると、正の勾配を持つギザギザの線が 200 単位の時間にわたってプロットされます。私が使用したプロットを使用するには、パッケージ「Rcell」が必要です。

壊滅的な状態を取り入れようとすると、私の困難が生じます。このコードに壊滅的な状態を組み込むにはどうすればよいですか? 私はこのようなものを想像しますが、構文をどのように翻訳するかわかりません:

step_generator <- function(col, rnum)

for (i in 2:length(col)

start in normal growth conditions (growth prob = A ; shrinkage prob = B)

if rnum[i] < C switch to catastrophic state (
              growth prob = A ; 
              shrinkage prob = E), 
else stay in normal growth conditions (i.e. if rnum[i] >= C)

stay in catastrophic state until rnum[i] < D, when this happens switch back to normal growth conditions. 

repeat through the entire 20000 rows

すべての助けに感謝します!

4

1 に答える 1

0

これは、乱数の1つの列で行われた私の質問に対する解決策です。複数の列に簡単に反復できます。

rnumbers <- data.frame(rnum = runif(20000, 0, 1)) #generates random #df 
dt <- c(.01) #time interval
A <- dt*1 #probability for growth
B <- dt*.25 # probability for shrinkage during growth phase
C <- dt*.03 # probability to enter catastrophe phase
D <- dt*.3 # probability to exit catastrophe phase and back into growth phase
E <- dt*5 # probability of shrinkage while in growth phase
rnumbers <- transform(rnumbers, cat = ifelse(rnum < C, .5, 0)) # generate column of times to enter catastrophe
rnumbers <- transform(rnumbers, rescue = ifelse(rnum < D, 1, 0)) # generate column of times to exit catastrophe
rnumbers <- transform(rnumbers, state = rowSums(rnumbers[,2:3])) # adds together two previous columns now (1.5 = catastrophe and 1 = growth and 0 = do nothing)
rnumbers[1,4] <- 1 #initialize in state 1 (growth phase)
rnumbers$state1 <- rnumbers$state # copies state into state 1
rnumbers$state1[rnumbers$state1 == 0] <- NA #makes any 0 in state1 NA
rnumbers$state1 <- na.locf(rnumbers$state1) # replaces NA with above row value
rnumbers <- transform(rnumbers, dynamics = ifelse(state1 < 1.5 , ifelse(rnum < B, -1, ifelse(rnum<A,1,0)), ifelse(rnum < A, 1, ifelse(rnum < E, -1,0))), interval = rep(dt)) # determines if polymer will grow or not
sum(rnumbers[,6])

rnumbers <- transform(rnumbers, time = cumsum(interval), step = cumsum(dynamics)) #cumulative sums of interval and growth or shrinkage
dev.new(width=5, height=5) #makes new plot
cplot(rnumbers,step ~ time, geom = "line") #plots length of polymer vs time

私以外の誰かが興味を持っている場合 - これは動的不安定性を考慮した微小管成長のシミュレーションです。ポリマーは最初の長さ (この場合は 0) から始まり、一定の確率で成長または収縮します。シミュレーション中の任意の時点で、特定の割合で大惨事が発生する可能性があります。これにより、収縮の確率が 20 倍になり、成長の確率が同じに保たれます。大惨事から抜け出すためには、一定の確率で起こる救助が必要です。

これは、大惨事が含まれていない場合のポリマーの成長の様子です。

library("Rcell")
require("zoo")
rnumbers <- data.frame(rnum = runif(20000, 0, 1))
dt <- c(.01)
A <- dt*1
B <- dt*.25
rnumbers <- transform(rnumbers, poly = ifelse(rnum < A, 1, ifelse(rnum < B, -1, 0)), dt = rep(dt))
rnumbers <- transform(rnumbers, time = cumsum(dt), step = cumsum(poly))
dev.new(width=5, height=5)
cplot(rnumbers,step ~ time, geom = "line")
于 2013-12-06T16:15:01.677 に答える