Hi Jim,
Thanks a lot for replying.
In fact I'm trying to run a simulation study that enables me to calculate
the Bayes risk of a sampling plan selected from progressively type-II
censored Weibull model. One of the steps involves evaluating the expected
test time, which is a rather complicated formula that involves nested
multiple summations where the counters of the summation signs are
dependent, that's why I thought of I should create the incomb() function
inside the loop, or may be I didn't figure out how to relate its arguments
to the ones inside the loop had I created it outside it. I'm trying to
create a matrix of all the possible combinations involved in the summations
and then use the apply() function on each row of that matrix. The problem
is that the code I wrote works perfectly well for rather small values of
the sample size,n, and the censoring number, m (for example, n=8,m=4),but
when n and m are increased (say, n=25,m=15) the code keeps on running for
days with no output. That's why I thought I should try to avoid explicit
loops as much as possible, so I did my best in this regard but still the
code takes too long to execute,(more than three days), thus, i believe
there must be something wrong.
Here's the full code:
library(pbapply)
f1 <- function(n, m) {
stopifnot(n > m)
r0 <- t(diff(combn(n-1, m-1)) - 1L)
r1 <- rep(seq(from=0, len=n-m+1), choose( seq(to=m-2, by=-1, len=n-m+1),
m-2))
cbind(r0[, ncol(r0):1, drop=FALSE], r1, deparse.level=0)
}
simpfun<- function (x,n,m,p,alpha,beta)
{
a<-factorial(n-m)/(prod((factorial(x)))*(factorial((n-m)- sum(x))))
b <- ((m-1):1)
c<- a*((p)^(sum(x)))*((1-p)^(((m-1)*(n-m))- sum(x%*%(as.matrix(b)))))
d <- n - cumsum(x) - (1:(m-1))
e<- n*(prod(d))*c
LD<-list()
for (i in 1:(m-1)) {
LD[[i]]<-seq(0,x[i],1)
}
LD[[m]]<-seq(0,(n-m-sum(x)),1)
LED<-expand.grid (LD)
LED<-as.matrix(LED)
store1<-numeric(nrow(LED))
for (j in 1:length(store1) )
{
incomb<-function(x,alpha,beta) {
g<-((-1)^(sum(LED[j,])))*(gamma((1/beta)+1))*((alpha)^(-(1/beta)))
h <- choose(x, LED[j,-m])
ik<-prod(h)*choose((n-m-sum(x)),LED[j,m])
lm<-cumsum(LED[j,-m]) + (1:(m-1))
plm<-prod(lm)
gil<-g*ik/(plm)
hlm<-numeric(sum(LED[j,])+(m-1))
dsa<-length(hlm)
for (i in 1:dsa)
{
ppp<- sum(LED[j,])+(m-1)
hlm[i]<-
(choose(ppp,i))*((-1)^(i))*((i+1)^((-1)*((1/beta)+1)))
}
shl<-gil*(sum(hlm)+1)
return (shl)
}
store1[j]<-incomb(x,alpha=0.2,beta=2)
}
val1<- sum(store1)*e
return(val1)
}
va<-pbapply(s,1,simpfun,n=6,m=4,p=0.3,alpha=0.2,beta=2)
EXP<-sum(va)
Any help would be greatly appreciated.
Thanks a lot for your time.
Best Regards,
Maram Salem
On 2 November 2015 at 00:27, jim holtman <jholt...@gmail.com> wrote:
> Why are you recreating the incomb function within the loop instead of
> defining it outside the loop? Also you are referencing several variables
> that are global (e.g., m & j); you should be passing these in as parameters
> to the function.
>
>
> Jim Holtman
> Data Munger Guru
>
> What is the problem that you are trying to solve?
> Tell me what you want to do, not how you want to do it.
>
> On Sun, Nov 1, 2015 at 7:31 AM, Maram SAlem <marammagdysa...@gmail.com>
> wrote:
>
>> Hi All,
>>
>> I'm writing a long code that takes long time to execute. So I used the
>> Rprof() function and found out that the function that takes about 80% of
>> the time is the incomb () fucntion (below), and this is most probably
>> because of the many explicit for() loops I'm using.
>>
>> n=18;m=4;p=0.3;alpha=0.2;beta=2
>> x=c(3,0,0)
>> LD<-list()
>> for (i in 1:(m-1)) {
>> LD[[i]]<-seq(0,x[i],1)
>> }
>> LD[[m]]<-seq(0,(n-m-sum(x)),1)
>> LED<-expand.grid (LD)
>> LED<-as.matrix(LED)
>> store1<-numeric(nrow(LED))
>> h<- numeric(m-1)
>> lm<- numeric(m-1)
>> for (j in 1:length(store1) )
>> {
>> incomb<-function(x,alpha,beta) {
>>
>> g<-((-1)^(sum(LED[j,])))*(gamma((1/beta)+1))*((alpha)^(-(1/beta)))
>> for (i in 1:(m-1)) {
>> h[i]<- choose(x[i],LED[j,i])
>> }
>> ik<-prod(h)*choose((n-m-sum(x)),LED[j,m])
>> for (i in 1:(m-1)) {
>> lm[i]<-(sum(LED[j,1:i])) + i
>> }
>> plm<-prod(lm)
>> gil<-g*ik/(plm)
>> hlm<-numeric(sum(LED[j,])+(m-1))
>> dsa<-length(hlm)
>> for (i in 1:dsa)
>> {
>> ppp<- sum(LED[j,])+(m-1)
>> hlm[i]<-
>> (choose(ppp,i))*((-1)^(i))*((i+1)^((-1)*((1/beta)+1)))
>> }
>> shl<-gil*(sum(hlm)+1)
>> return (shl)
>> }
>> store1[j]<-incomb(x,alpha=0.2,beta=2)
>> }
>>
>>
>> I'm trying to use alternatives (for ex. to vectorize things) to the
>> explicit for() loops, but things don't work out.
>>
>> Any suggestions that can help me to speed up the execution of the incomb()
>> function are much appreciated.
>>
>> Thanks a lot in advance.
>>
>> Maram Salem
>>
>> [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
