Uwe Ligges
2005-Mar-02 17:03 UTC
[R] Warning: number of items to replace is not a multiple of replacement length
Lyle W. Konigsberg wrote:> I feel like a complete dolt, as I know this question has been asked by > others on a fairly regular basis, but I'm going in circles trying to get > the following to work: > > id.prob<-function (tt) > { > library(mvtnorm) > #===========================> Makeham<-function(tt) > { > a2=0.030386513 > a3=0.006688287 > b3=0.039047537 > t<-tt-20 > h.t<-a2+a3*exp(b3*t) > S.t<-exp(-a2*t+a3/b3*(1-exp(b3*t))) > return(S.t*h.t) > } > #==========================> trans<-function (age) > { > indic=c(2,2) > lage=log(age) > xpars=c(2.1862908,7.5528077,8.5806697,2.3319461,8.8959507,9.1380187,0.3720293) > > beta<-c(0,0) > alpha1<-c(0,0) > alpha2<-c(0,0) > beta[1]<-xpars[1] > alpha1[1]<-xpars[2] > alpha2[1]<-xpars[3] > beta[2]<-xpars[4] > alpha1[2]<-xpars[5] > alpha2[2]<-xpars[6] > corr<-matrix(rep(0,4),nc=2) > corr[1,1]=1 > corr[2,2]=1 > corr[1,2]=xpars[7] > corr[2,1]=corr[1,2] > > LLK<-0 > L<-c(0,0) > R<-c(0,0) > > for(j in 1:2){ > if(indic[j]==1){ > L[j]<--Inf > R[j]<-alpha1[j]-beta[j]*lage} > if(indic[j]==2){ > L[j]<-alpha1[j]-beta[j]*lage > R[j]<-alpha2[j]-beta[j]*lage} > if(indic[j]==3){ > L[j]<-alpha2[j]-beta[j]*lage > R[j]<-Inf} > } > > LLK<-pmvnorm(lower=L,upper=R,corr=corr)[1] > > return(LLK) > } > > #==========================> > tot<-Makeham(tt)*trans(tt) > return(tot) > } > > > And then: > > > id.prob(20) > [1] 0.0001417763 > > id.prob(21) > [1] 0.00018078 > > id.prob(20:21) > [1] 0.0001417763 0.0001375794 > Warning messages: > 1: number of items to replace is not a multiple of replacement length > 2: number of items to replace is not a multiple of replacement length > 3: number of items to replace is not a multiple of replacement length > 4: number of items to replace is not a multiple of replacement lengthPlease read some introductory literature regarding the S / R language! Your code is not vectorizable at all. In tot <- Makeham(tt)*trans(tt) Makeham has length 2 and works (I guess), but in trans(tt) you have lage=log(age) which has length 2, but you get also length 2 in R[j]<-alpha1[j]-beta[j]*lage} but assign it to a scalar value.> Ultimately I need to be able to integrate this function, which isn't > going to happen (correctly) if I've messed-up writing the function. > Thanks for any advice on this. > > Regards, > Lyle W. Konigsberg > Professor of Anthropology > University of Tennessee > Knoxville, TN > LyleK at utk.edu > > ______________________________________________ > R-help at stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html
Lyle W. Konigsberg
2005-Mar-02 18:58 UTC
[R] Warning: number of items to replace is not a multiple of replacement length
I feel like a complete dolt, as I know this question has been asked by
others on a fairly regular basis, but I'm going in circles trying to get
the following to work:
id.prob<-function (tt)
{
library(mvtnorm)
#===========================Makeham<-function(tt)
{
a2=0.030386513
a3=0.006688287
b3=0.039047537
t<-tt-20
h.t<-a2+a3*exp(b3*t)
S.t<-exp(-a2*t+a3/b3*(1-exp(b3*t)))
return(S.t*h.t)
}
#==========================trans<-function (age)
{
indic=c(2,2)
lage=log(age)
xpars=c(2.1862908,7.5528077,8.5806697,2.3319461,8.8959507,9.1380187,0.3720293)
beta<-c(0,0)
alpha1<-c(0,0)
alpha2<-c(0,0)
beta[1]<-xpars[1]
alpha1[1]<-xpars[2]
alpha2[1]<-xpars[3]
beta[2]<-xpars[4]
alpha1[2]<-xpars[5]
alpha2[2]<-xpars[6]
corr<-matrix(rep(0,4),nc=2)
corr[1,1]=1
corr[2,2]=1
corr[1,2]=xpars[7]
corr[2,1]=corr[1,2]
LLK<-0
L<-c(0,0)
R<-c(0,0)
for(j in 1:2){
if(indic[j]==1){
L[j]<--Inf
R[j]<-alpha1[j]-beta[j]*lage}
if(indic[j]==2){
L[j]<-alpha1[j]-beta[j]*lage
R[j]<-alpha2[j]-beta[j]*lage}
if(indic[j]==3){
L[j]<-alpha2[j]-beta[j]*lage
R[j]<-Inf}
}
LLK<-pmvnorm(lower=L,upper=R,corr=corr)[1]
return(LLK)
}
#==========================
tot<-Makeham(tt)*trans(tt)
return(tot)
}
And then:
> id.prob(20)
[1] 0.0001417763
> id.prob(21)
[1] 0.00018078
> id.prob(20:21)
[1] 0.0001417763 0.0001375794
Warning messages:
1: number of items to replace is not a multiple of replacement length
2: number of items to replace is not a multiple of replacement length
3: number of items to replace is not a multiple of replacement length
4: number of items to replace is not a multiple of replacement length
Ultimately I need to be able to integrate this function, which isn't going
to happen (correctly) if I've messed-up writing the function. Thanks for
any advice on this.
Regards,
Lyle W. Konigsberg
Professor of Anthropology
University of Tennessee
Knoxville, TN
LyleK at utk.edu