From 17602bc2df5ab65aea88eb90bf07e575f1e49e07 Mon Sep 17 00:00:00 2001 From: "N. Brouard" Date: Fri, 13 Jun 2003 07:45:28 +0000 Subject: [PATCH] We decided to release, in June 2003, a version 0.8a1 (May 2002) of former version 0.8a to correct a bug highlighted by Chris Jackson. It can be useful until version 0.9x will be released. In version 0.8a, likelihood was computed as if date of death was unknown. Death was treated as any other health state, i.e the date of the interview describes the actual state and not the date of a change in health state. It means that in the early version of IMaCh we were considering that at each interview the current state was recorded (healthy, disable or death). And IMaCh was correct. But, when we introduced the exact date of death then we should have modified the contribution of an exact (more precise) date of death to the likelihood. This new contribution is smaller and very dependent of the step unit stepm. It is no more the probability to die between last interview and month of death (which is the way version 0.8a of IMaCh was doing) but the probability to survive from last interview up to one month before death multiplied by the probability to die within a month. Thanks to Chris Jackson for correcting this bug. Former versions increased mortality artificially. The bad side is that we add another loop which slows down the processing. The difference can be up to 10% lower mortality. --- src/imach.c | 42 ++++++++++++++++++++++++++++++++++++++---- 1 file changed, 38 insertions(+), 4 deletions(-) diff --git a/src/imach.c b/src/imach.c index 636f597..3efca52 100644 --- a/src/imach.c +++ b/src/imach.c @@ -870,6 +870,7 @@ double func( double *x) double **out; double sw; /* Sum of weights */ double lli; /* Individual log likelihood */ + int s1, s2; long ipmx; /*extern weight */ /* We are differentiating ll according to initial status */ @@ -884,7 +885,10 @@ double func( double *x) for (k=1; k<=cptcovn;k++) cov[2+k]=covar[Tvar[k]][i]; for(mi=1; mi<= wav[i]-1; mi++){ for (ii=1;ii<=nlstate+ndeath;ii++) - for (j=1;j<=nlstate+ndeath;j++) oldm[ii][j]=(ii==j ? 1.0 : 0.0); + for (j=1;j<=nlstate+ndeath;j++){ + oldm[ii][j]=(ii==j ? 1.0 : 0.0); + savm[ii][j]=(ii==j ? 1.0 : 0.0); + } for(d=0; d nlstate){ + /* i.e. if s2 is a death state and if the date of death is known then the contribution + to the likelihood is the probability to die between last step unit time and current + step unit time, which is also the differences between probability to die before dh + and probability to die before dh-stepm . + In version up to 0.92 likelihood was computed + as if date of death was unknown. Death was treated as any other + health state: the date of the interview describes the actual state + and not the date of a change in health state. The former idea was + to consider that at each interview the state was recorded + (healthy, disable or death) and IMaCh was corrected; but when we + introduced the exact date of death then we should have modified + the contribution of an exact death to the likelihood. This new + contribution is smaller and very dependent of the step unit + stepm. It is no more the probability to die between last interview + and month of death but the probability to survive from last + interview up to one month before death multiplied by the + probability to die within a month. Thanks to Chris + Jackson for correcting this bug. Former versions increased + mortality artificially. The bad side is that we add another loop + which slows down the processing. The difference can be up to 10% + lower mortality. + */ + lli=log(out[s1][s2] - savm[s1][s2]); + }else{ + lli=log(out[s1][s2]); /* or lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]); */ + /* printf(" %f ",out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/ + } ipmx +=1; sw += weight[i]; ll[s[mw[mi][i]][i]] += 2*weight[i]*lli; + /*printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d lli=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],lli,weight[i],out[s1][s2],savm[s1][s2]);*/ } /* end of wave */ } /* end of individual */ for(k=1,l=0.; k<=nlstate; k++) l += ll[k]; /* printf("l1=%f l2=%f ",ll[1],ll[2]); */ l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */ + /*exit(0);*/ return -l; } @@ -2645,7 +2679,7 @@ int main(int argc, char *argv[]) double dateprev1, dateprev2,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2; - char version[80]="Imach version 0.8a, May 2002, INED-EUROREVES "; + char version[80]="Imach version 0.8a1, June 2003, INED-EUROREVES "; char *alph[]={"a","a","b","c","d","e"}, str[4]; -- 2.43.0