-1

我有以下等式，我必须用 R 来写它。在此处输入图像描述

我在使用这种编码时面临的主要挑战是，对于每个 i 和每个 k，我需要一个组合总和（分子中的二项式）。这里对于某些 i 和 k，min{(nk),Mi} 可以为 0，特别是当 k=n 时。最重要的是，循环不能从 0 开始，但我需要它！

出于您的考虑，这是我的代码和数据（d1）。你会看到我需要在 line 之后从 0 开始循环sum<-0，这是我的主要问题。你能看到问题并更正代码吗？需要从 0 开始的循环怎么办？

n<-4

n^2

id<-1:16

r<-c(1,1,1,1,0,1,1,1,0,1,0,1,1,1,1,1)

tr<-rep(0,n^2)
for(i in 1:n^2){
tr[i]<-ifelse(r[i]==1,rexp(1,1/1),0)
}

t0<-rep(0,n^2)
for(i in 1:n^2){
t0[i]<-ifelse(r[i]==0,rexp(1,1/1.5),0)
}

#Total number of subjects who cannot get B1:

M<-sum(r==0)    #If there were cenoring then M<-sum(r==0 & tr>cenc)

d<-data.frame(id,r,t0,tr)

d1<-d[tr>0,]
d1

d1<-d1[order(d1$tr),]
d1

d1$rank<-1:length(d1$tr)
d1


###Calculating the probability of getting B1 for each subject with r[i]=1:

d1$prob<-rep(0,length(d1$tr))
for(i in 1:length(d1$tr)){  #loop i begins
Mi<-sum(d1$tr[i]>t0[t0>0])
for(k in 1:n){  #loop k begins
sum<-0
for(m in 0:min(n-k,Mi)){    #loop m begins
sum<-sum+choose(Mi,m)*choose(n^2-i-Mi,n-k-m)
}   #loop m ends.
d1$prob[i]<-d1$prob[i]+choose(i-1,k-1)*sum/choose(n^2,n)
}   #loop k ends.
}   #loop i ends.
d1$prob<-d1$prob*1/n
4

1 回答 1

2

我看不出您的原始代码有什么问题？您不能索引从零开始的向量，但它与循环中迭代器的范围无关。无论如何，这应该给你相同的答案，但避免了最里面的循环：

d1$prob<-rep(0,length(d1$tr))
for(i in 1:length(d1$tr)){  #loop i begins
  Mi<-sum(d1$tr[i]>t0[t0>0])
  for(k in 1:n){  #loop k begins
    # function choose is vectorized, so you can compute all 
    # binomial terms at once given k
    SUM<-sum(choose(Mi,0:min(n-k,Mi))*choose(n^2-i-Mi,n-k-(0:min(n-k,Mi))))
    d1$prob[i]<-d1$prob[i]+choose(i-1,k-1)*SUM
  }   #loop k ends.
}   #loop i ends.
d1$prob<-(d1$prob/choose(n^2,n))*1/n 
#as choose(n^2,n) does not depend on i nor k, 
#you can make the division after the loops for all elements at once
于 2013-03-24T19:30:44.883 回答