Marion Keast
2015-Sep-23 09:52 UTC
[R] Issues fitting Gompertz-Makeham model to mortality data
Good evening, I have a very basic knowledge of how to use R and am struggling with my code to fit a Gompertz-Makeham model to data. The code and error message is as follows:> data <- hmd.mx("AUS","mkeast93 at hotmail.com","1Mrkqazwsx", "country")Warning message: In hmd.mx("AUS", "mkeast93 at hotmail.com", "1Mrkqazwsx", "country") : NAs introduced by coercion> summary(data)Mortality data for country Series: female male total Years: 1921 - 2011 Ages: 0 - 110> m.qx <- data$rate$male[1:111,1:90] > f.qx <-data$rate$female[1:111,1:90] > age<-data$age[21:90] > > n=111 > m=90 > beta1=(1:m)*0 > beta2=(1:m)*0 > beta3=(1:m)*0 > qhat.mx=array(0,dim=c(n,m)) > qhat.fx=array(0,dim=c(n,m)) > > #j is ages 0 - 110 > #i is years 1921 - 2011 > > #MALES > for (i in 1:m)+ { + model1 <-fitGM(,m.qx[1:n,i]) + m1 <- model1[1] + m2 <- model1[2] + m3 <- model1[3] + beta1[i]=m1 + beta2[i]=m2 + beta3[i]=m3 + qhat.mx[,i]=GompertzMakeham(beta1[i],beta2[i],beta3[i],age) + } Error in qhat.mx[, i] = GompertzMakeham(beta1[i], beta2[i], beta3[i], : number of items to replace is not a multiple of replacement length> > MAPE.m=sum(abs(qhat.mx-m.qx)/m.qx)/(n*m) > MaPE.m=sum((qhat.mx-m.qx)/m.qx)/(n*m)Please advise in regards to fixing errors. Furthermore, if I would like to fit MAPE to each individual age instead of the population, is there a loop I can use to achieve this. Thank you very much for your time and assistance. Marion [[alternative HTML version deleted]]