/* $Id$
$State$
$Log$
+ Revision 1.213 2015/12/11 18:22:17 brouard
+ Summary: 0.98r4
+
Revision 1.212 2015/11/21 12:47:24 brouard
Summary: minor typo
double **oldm, **newm, **savm; /* Working pointers to matrices */
double **oldms, **newms, **savms; /* Fixed working pointers to matrices */
/*FILE *fic ; */ /* Used in readdata only */
-FILE *ficpar, *ficparo,*ficres, *ficresp, *ficresphtm, *ficrespl, *ficrespij, *ficrest,*ficresf,*ficrespop;
+FILE *ficpar, *ficparo,*ficres, *ficresp, *ficresphtm, *ficresphtmfr, *ficrespl, *ficrespij, *ficrest,*ficresf,*ficrespop;
FILE *ficlog, *ficrespow;
int globpr=0; /* Global variable for printing or not */
double fretone; /* Only one call to likelihood */
double **out, cov[NCOVMAX+1];
double **newm;
double agexact;
+ double agebegin, ageend;
/* Hstepm could be zero and should return the unit matrix */
for (i=1;i<=nlstate+ndeath;i++)
newm=savm;
/* Covariates have to be included here again */
cov[1]=1.;
- agexact=age+((h-1)*hstepm + (d-1))*stepm/YEARM;
+ agexact=age+((h-1)*hstepm + (d-1))*stepm/YEARM; /* age just before transition */
cov[2]=agexact;
if(nagesqr==1)
cov[3]= agexact*agexact;
int s1, s2;
double bbh, survp;
double agexact;
+ double agebegin, ageend;
/*extern weight */
/* We are differentiating ll according to initial status */
/* for (i=1;i<=npar;i++) printf("%f ", x[i]);*/
oldm[ii][j]=(ii==j ? 1.0 : 0.0);
savm[ii][j]=(ii==j ? 1.0 : 0.0);
}
- for(d=0; d<dh[mi][i]; d++){
+
+ agebegin=agev[mw[mi][i]][i]; /* Age at beginning of effective wave */
+ ageend=agev[mw[mi][i]][i] + (dh[mi][i])*stepm/YEARM; /* Age at end of effective wave and at the end of transition */
+ for(d=0; d<dh[mi][i]; d++){ /* Delay between two effective waves */
+ /*dh[m][i] or dh[mw[mi][i]][i] is the delay between two effective waves m=mw[mi][i]
+ and mw[mi+1][i]. dh depends on stepm.*/
newm=savm;
agexact=agev[mw[mi][i]][i]+d*stepm/YEARM;
cov[2]=agexact;
ll[s[mw[mi][i]][i]] += 2*weight[i]*lli;
/*printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],out[s1][s2],savm[s1][s2]); */
if(globpr){
- fprintf(ficresilk,"%9ld %6.1f %6d %2d %2d %2d %2d %3d %11.6f %8.4f %8.3f\
+ fprintf(ficresilk,"%9ld %6.1f %6.1f %6d %2d %2d %2d %2d %3d %11.6f %8.4f %8.3f\
%11.6f %11.6f %11.6f ", \
- num[i], agexact, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw,
+ num[i], agebegin, ageend, i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],weight[i]*gipmx/gsw,
2*weight[i]*lli,out[s1][s2],savm[s1][s2]);
for(k=1,llt=0.,l=0.; k<=nlstate; k++){
llt +=ll[k]*gipmx/gsw;
printf("Problem with resultfile: %s\n", fileresilk);
fprintf(ficlog,"Problem with resultfile: %s\n", fileresilk);
}
- fprintf(ficresilk, "#individual(line's_record) count age s1 s2 wave# effective_wave# number_of_matrices_product pij weight weight/gpw -2ln(pij)*weight 0pij_x 0pij_(x-stepm) cumulating_loglikeli_by_health_state(reweighted=-2ll*weightXnumber_of_contribs/sum_of_weights) and_total\n");
- fprintf(ficresilk, "#num_i age i s1 s2 mi mw dh likeli weight %%weight 2wlli out sav ");
+ fprintf(ficresilk, "#individual(line's_record) count ageb ageend s1 s2 wave# effective_wave# number_of_matrices_product pij weight weight/gpw -2ln(pij)*weight 0pij_x 0pij_(x-stepm) cumulating_loglikeli_by_health_state(reweighted=-2ll*weightXnumber_of_contribs/sum_of_weights) and_total\n");
+ fprintf(ficresilk, "#num_i ageb agend i s1 s2 mi mw dh likeli weight %%weight 2wlli out sav ");
/* i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],2*weight[i]*lli,out[s1][s2],savm[s1][s2]); */
for(k=1; k<=nlstate; k++)
fprintf(ficresilk," -2*gipw/gsw*weight*ll[%d]++",k);
}
/************ Frequencies ********************/
-void freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, int *Tvaraff, int **nbcode, int *ncodemax,double **mint,double **anint, char strstart[])
+void freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, \
+ int *Tvaraff, int **nbcode, int *ncodemax,double **mint,double **anint, char strstart[],\
+ int firstpass, int lastpass, int stepm, int weightopt, char model[])
{ /* Some frequencies */
int i, m, jk, j1, bool, z1,j;
+ int mi; /* Effective wave */
int first;
double ***freq; /* Frequencies */
double *pp, **prop;
double pos,posprop, k2, dateintsum=0,k2cpt=0;
- char fileresp[FILENAMELENGTH], fileresphtm[FILENAMELENGTH];
-
+ char fileresp[FILENAMELENGTH], fileresphtm[FILENAMELENGTH], fileresphtmfr[FILENAMELENGTH];
+ double agebegin, ageend;
+
pp=vector(1,nlstate);
prop=matrix(1,nlstate,iagemin,iagemax+3);
strcpy(fileresp,"P_");
fprintf(ficlog,"Problem with prevalence resultfile: %s\n", fileresp);
exit(0);
}
- printf("Problem with prevalence resultfile: %s\n", fileresp);
+
strcpy(fileresphtm,subdirfext(optionfilefiname,"PHTM_",".htm"));
if((ficresphtm=fopen(fileresphtm,"w"))==NULL) {
printf("Problem with prevalence HTM resultfile '%s' with errno='%s'\n",fileresphtm,strerror(errno));
fflush(ficlog);
exit(70);
}
+ else{
+ fprintf(ficresphtm,"<html><head>\n<title>IMaCh PHTM_ %s</title></head>\n <body><font size=\"2\">%s <br> %s</font> \
+<hr size=\"2\" color=\"#EC5E5E\"> \n\
+Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s<br>\n",\
+ fileresphtm,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model);
+ }
+ fprintf(ficresphtm,"Current page is file <a href=\"%s\">%s</a><br>\n\n<h4>Frequencies and prevalence by age at begin of transition</h4>\n",fileresphtm, fileresphtm);
+
+ strcpy(fileresphtmfr,subdirfext(optionfilefiname,"PHTMFR_",".htm"));
+ if((ficresphtmfr=fopen(fileresphtmfr,"w"))==NULL) {
+ printf("Problem with frequency table HTM resultfile '%s' with errno='%s'\n",fileresphtmfr,strerror(errno));
+ fprintf(ficlog,"Problem with frequency table HTM resultfile '%s' with errno='%s'\n",fileresphtmfr,strerror(errno));
+ fflush(ficlog);
+ exit(70);
+ }
+ else{
+ fprintf(ficresphtmfr,"<html><head>\n<title>IMaCh PHTM_Frequency table %s</title></head>\n <body><font size=\"2\">%s <br> %s</font> \
+<hr size=\"2\" color=\"#EC5E5E\"> \n\
+Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=1+age+%s<br>\n",\
+ fileresphtmfr,version,fullversion,title,datafile,firstpass,lastpass,stepm, weightopt, model);
+ }
+ fprintf(ficresphtmfr,"Current page is file <a href=\"%s\">%s</a><br>\n\n<h4>Frequencies of all effective transitions by age at begin of transition </h4>Unknown status is -1<br/>\n",fileresphtmfr, fileresphtmfr);
+
freq= ma3x(-5,nlstate+ndeath,-5,nlstate+ndeath,iagemin,iagemax+3);
j1=0;
first=1;
- /* for(k1=1; k1<=j ; k1++){ */ /* Loop on covariates */
- /* for(i1=1; i1<=ncodemax[k1];i1++){ */ /* Now it is 2 */
- /* j1++; */
- for (j1 = 1; j1 <= (int) pow(2,cptcoveff); j1++){
+ for (j1 = 1; j1 <= (int) pow(2,cptcoveff); j1++){ /* Loop on covariates combination */
/*printf("cptcoveff=%d Tvaraff=%d", cptcoveff,Tvaraff[1]);
scanf("%d", i);*/
for (i=-5; i<=nlstate+ndeath; i++)
dateintsum=0;
k2cpt=0;
- for (i=1; i<=imx; i++) {
+ for (i=1; i<=imx; i++) { /* For each individual i */
bool=1;
if (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
for (z1=1; z1<=cptcoveff; z1++)
/* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/
}
} /* cptcovn > 0 */
-
+
if (bool==1){
- for(m=firstpass; m<=lastpass; m++){
- k2=anint[m][i]+(mint[m][i]/12.);
- /*if ((k2>=dateprev1) && (k2<=dateprev2)) {*/
- if(agev[m][i]==0) agev[m][i]=iagemax+1;
- if(agev[m][i]==1) agev[m][i]=iagemax+2;
- if (s[m][i]>0 && s[m][i]<=nlstate)
- prop[s[m][i]][(int)agev[m][i]] += weight[i];
+ /* for(m=firstpass; m<=lastpass; m++){ */
+ for(mi=1; mi<wav[i];mi++){
+ m=mw[mi][i];
+ /* dh[m][i] or dh[mw[mi][i]][i] is the delay between two effective (mi) waves m=mw[mi][i]
+ and mw[mi+1][i]. dh depends on stepm. */
+ agebegin=agev[m][i]; /* Age at beginning of wave before transition*/
+ ageend=agev[m][i]+(dh[m][i])*stepm/YEARM; /* Age at end of wave and transition */
+ if(m >=firstpass && m <=lastpass){
+ k2=anint[m][i]+(mint[m][i]/12.);
+ /*if ((k2>=dateprev1) && (k2<=dateprev2)) {*/
+ if(agev[m][i]==0) agev[m][i]=iagemax+1; /* All ages equal to 0 are in iagemax+1 */
+ if(agev[m][i]==1) agev[m][i]=iagemax+2; /* All ages equal to 1 are in iagemax+2 */
+ if (s[m][i]>0 && s[m][i]<=nlstate) /* If status at wave m is known and a live state */
+ prop[s[m][i]][(int)agev[m][i]] += weight[i]; /* At age of beginning of transition, where status is known */
if (m<lastpass) {
- freq[s[m][i]][s[m+1][i]][(int)agev[m][i]] += weight[i];
- freq[s[m][i]][s[m+1][i]][iagemax+3] += weight[i];
+ /* if(s[m][i]==4 && s[m+1][i]==4) */
+ /* printf(" num=%ld m=%d, i=%d s1=%d s2=%d agev at m=%d\n", num[i], m, i,s[m][i],s[m+1][i], (int)agev[m][i]); */
+ if(s[m][i]==-1)
+ printf(" num=%ld m=%d, i=%d s1=%d s2=%d agev at m=%d agebegin=%.2f ageend=%.2f, agemed=%d\n", num[i], m, i,s[m][i],s[m+1][i], (int)agev[m][i],agebegin, ageend, (int)((agebegin+ageend)/2.));
+ freq[s[m][i]][s[m+1][i]][(int)agev[m][i]] += weight[i]; /* At age of beginning of transition, where status is known */
+ /* freq[s[m][i]][s[m+1][i]][(int)((agebegin+ageend)/2.)] += weight[i]; */
+ freq[s[m][i]][s[m+1][i]][iagemax+3] += weight[i]; /* Total is in iagemax+3 *//* At age of beginning of transition, where status is known */
}
-
- if ((agev[m][i]>1) && (agev[m][i]< (iagemax+3)) && (anint[m][i]!=9999) && (mint[m][i]!=99)) {
- dateintsum=dateintsum+k2;
- k2cpt++;
- /* printf("i=%ld dateintmean = %lf dateintsum=%lf k2cpt=%lf k2=%lf\n",i, dateintsum/k2cpt, dateintsum,k2cpt, k2); */
- }
- /*}*/
+ }
+ if ((agev[m][i]>1) && (agev[m][i]< (iagemax+3)) && (anint[m][i]!=9999) && (mint[m][i]!=99)) {
+ dateintsum=dateintsum+k2;
+ k2cpt++;
+ /* printf("i=%ld dateintmean = %lf dateintsum=%lf k2cpt=%lf k2=%lf\n",i, dateintsum/k2cpt, dateintsum,k2cpt, k2); */
+ }
+ /*}*/
} /* end m */
} /* end bool */
} /* end i = 1 to imx */
pstamp(ficresp);
if (cptcovn>0) {
fprintf(ficresp, "\n#********** Variable ");
- fprintf(ficresphtm, "\n<h3>********** Variable ");
+ fprintf(ficresphtm, "\n<br/><br/><h3>********** Variable ");
+ fprintf(ficresphtmfr, "\n<br/><br/><h3>********** Variable ");
for (z1=1; z1<=cptcoveff; z1++){
fprintf(ficresp, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
fprintf(ficresphtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
+ fprintf(ficresphtmfr, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
}
fprintf(ficresp, "**********\n#");
- fprintf(ficresphtm, "**********</h3>\n#");
+ fprintf(ficresphtm, "**********</h3>\n");
+ fprintf(ficresphtmfr, "**********</h3>\n");
fprintf(ficlog, "\n#********** Variable ");
for (z1=1; z1<=cptcoveff; z1++) fprintf(ficlog, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]);
- fprintf(ficlog, "**********\n#");
+ fprintf(ficlog, "**********\n");
}
- fprintf(ficresphtm,"<table style=\"text-align:center; border: 1px solid\"><th></th>");
+ fprintf(ficresphtm,"<table style=\"text-align:center; border: 1px solid\">");
for(i=1; i<=nlstate;i++) {
fprintf(ficresp, " Age Prev(%d) N(%d) N",i,i);
fprintf(ficresphtm, "<th>Age</th><th>Prev(%d)</th><th>N(%d)</th><th>N</th>",i,i);
fprintf(ficresp, "\n");
fprintf(ficresphtm, "\n");
+ /* Header of frequency table by age */
+ fprintf(ficresphtmfr,"<table style=\"text-align:center; border: 1px solid\">");
+ fprintf(ficresphtmfr,"<th>Age</th> ");
+ for(jk=-1; jk <=nlstate+ndeath; jk++){
+ for(m=-1; m <=nlstate+ndeath; m++){
+ if(jk!=0 && m!=0)
+ fprintf(ficresphtmfr,"<th>%d%d</th> ",jk,m);
+ }
+ }
+ fprintf(ficresphtmfr, "\n");
+
+ /* For each age */
for(i=iagemin; i <= iagemax+3; i++){
fprintf(ficresphtm,"<tr>");
- if(i==iagemax+3){
+ if(i==iagemax+1){
+ fprintf(ficlog,"1");
+ fprintf(ficresphtmfr,"<tr><th>0</th> ");
+ }else if(i==iagemax+2){
+ fprintf(ficlog,"0");
+ fprintf(ficresphtmfr,"<tr><th>Unknown</th> ");
+ }else if(i==iagemax+3){
fprintf(ficlog,"Total");
+ fprintf(ficresphtmfr,"<tr><th>Total</th> ");
}else{
if(first==1){
first=0;
printf("See log file for details...\n");
}
+ fprintf(ficresphtmfr,"<tr><th>%d</th> ",i);
fprintf(ficlog,"Age %d", i);
}
for(jk=1; jk <=nlstate ; jk++){
}
}
- for(jk=-1; jk <=nlstate+ndeath; jk++)
- for(m=-1; m <=nlstate+ndeath; m++)
- if(freq[jk][m][i] !=0 ) {
- if(first==1)
- printf(" %d%d=%.0f",jk,m,freq[jk][m][i]);
+ for(jk=-1; jk <=nlstate+ndeath; jk++){
+ for(m=-1; m <=nlstate+ndeath; m++){
+ if(freq[jk][m][i] !=0 ) { /* minimizing output */
+ if(first==1){
+ printf(" %d%d=%.0f",jk,m,freq[jk][m][i]);
+ }
fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][i]);
}
+ if(jk!=0 && m!=0)
+ fprintf(ficresphtmfr,"<td>%.0f</td> ",freq[jk][m][i]);
+ }
+ }
+ fprintf(ficresphtmfr,"</tr>\n ");
if(i <= iagemax){
fprintf(ficresp,"\n");
fprintf(ficresphtm,"</tr>\n");
fprintf(ficlog,"\n");
} /* end loop i */
fprintf(ficresphtm,"</table>\n");
+ fprintf(ficresphtmfr,"</table>\n");
/*}*/
} /* end j1 */
dateintmean=dateintsum/k2cpt;
fclose(ficresp);
fclose(ficresphtm);
+ fclose(ficresphtmfr);
free_ma3x(freq,-5,nlstate+ndeath,-5,nlstate+ndeath, iagemin, iagemax+3);
free_vector(pp,1,nlstate);
free_matrix(prop,1,nlstate,iagemin, iagemax+3);
*/
int i, m, jk, j1, bool, z1,j;
+ int mi; /* Effective wave */
+ int iage;
+ double agebegin, ageend;
double **prop;
double posprop;
first=1;
for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){
- /*for(i1=1; i1<=ncodemax[k1];i1++){
- j1++;*/
-
- for (i=1; i<=nlstate; i++)
- for(m=iagemin; m <= iagemax+3; m++)
- prop[i][m]=0.0;
-
- for (i=1; i<=imx; i++) { /* Each individual */
- bool=1;
- if (cptcovn>0) {
- for (z1=1; z1<=cptcoveff; z1++)
- if (covar[Tvaraff[z1]][i]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)])
- bool=0;
- }
- if (bool==1) {
- for(m=firstpass; m<=lastpass; m++){/* Other selection (we can limit to certain interviews*/
+ for (i=1; i<=nlstate; i++)
+ for(iage=iagemin; iage <= iagemax+3; iage++)
+ prop[i][iage]=0.0;
+
+ for (i=1; i<=imx; i++) { /* Each individual */
+ bool=1;
+ if (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
+ for (z1=1; z1<=cptcoveff; z1++)
+ if (covar[Tvaraff[z1]][i]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)])
+ bool=0;
+ }
+ if (bool==1) {
+ /* for(m=firstpass; m<=lastpass; m++){/\* Other selection (we can limit to certain interviews*\/ */
+ for(mi=1; mi<wav[i];mi++){
+ m=mw[mi][i];
+ agebegin=agev[m][i]; /* Age at beginning of wave before transition*/
+ /* ageend=agev[m][i]+(dh[m][i])*stepm/YEARM; /\* Age at end of wave and transition *\/ */
+ if(m >=firstpass && m <=lastpass){
y2=anint[m][i]+(mint[m][i]/12.); /* Fractional date in year */
if ((y2>=dateprev1) && (y2<=dateprev2)) { /* Here is the main selection (fractional years) */
if(agev[m][i]==0) agev[m][i]=iagemax+1;
if(agev[m][i]==1) agev[m][i]=iagemax+2;
if((int)agev[m][i] <iagemin || (int)agev[m][i] >iagemax+3) printf("Error on individual =%d agev[m][i]=%f m=%d\n",i, agev[m][i],m);
- if (s[m][i]>0 && s[m][i]<=nlstate) {
+ if (s[m][i]>0 && s[m][i]<=nlstate) {
/*if(i>4620) printf(" i=%d m=%d s[m][i]=%d (int)agev[m][i]=%d weight[i]=%f prop=%f\n",i,m,s[m][i],(int)agev[m][m],weight[i],prop[s[m][i]][(int)agev[m][i]]);*/
- prop[s[m][i]][(int)agev[m][i]] += weight[i];
- prop[s[m][i]][iagemax+3] += weight[i];
- }
- }
+ prop[s[m][i]][(int)agev[m][i]] += weight[i];/* At age of beginning of transition, where status is known */
+ prop[s[m][i]][iagemax+3] += weight[i];
+ } /* end valid statuses */
+ } /* end selection of dates */
} /* end selection of waves */
- }
- }
- for(i=iagemin; i <= iagemax+3; i++){
- for(jk=1,posprop=0; jk <=nlstate ; jk++) {
- posprop += prop[jk][i];
- }
-
- for(jk=1; jk <=nlstate ; jk++){
- if( i <= iagemax){
- if(posprop>=1.e-5){
- probs[i][jk][j1]= prop[jk][i]/posprop;
- } else{
- if(first==1){
- first=0;
- printf("Warning Observed prevalence probs[%d][%d][%d]=%lf because of lack of cases\nSee others on log file...\n",jk,i,j1,probs[i][jk][j1]);
- }
+ } /* end effective waves */
+ } /* end bool */
+ }
+ for(i=iagemin; i <= iagemax+3; i++){
+ for(jk=1,posprop=0; jk <=nlstate ; jk++) {
+ posprop += prop[jk][i];
+ }
+
+ for(jk=1; jk <=nlstate ; jk++){
+ if( i <= iagemax){
+ if(posprop>=1.e-5){
+ probs[i][jk][j1]= prop[jk][i]/posprop;
+ } else{
+ if(first==1){
+ first=0;
+ printf("Warning Observed prevalence probs[%d][%d][%d]=%lf because of lack of cases\nSee others on log file...\n",jk,i,j1,probs[i][jk][j1]);
}
- }
- }/* end jk */
- }/* end i */
+ }
+ }
+ }/* end jk */
+ }/* end i */
/*} *//* end i1 */
} /* end j1 */
int i, mi, m;
/* int j, k=0,jk, ju, jl,jmin=1e+5, jmax=-1;
double sum=0., jmean=0.;*/
- int first;
+ int first, firstwo;
int j, k=0,jk, ju, jl;
double sum=0.;
first=0;
+ firstwo=0;
jmin=100000;
jmax=-1;
jmean=0.;
- for(i=1; i<=imx; i++){
+ for(i=1; i<=imx; i++){ /* For simple cases and if state is death */
mi=0;
m=firstpass;
- while(s[m][i] <= nlstate){
+ while(s[m][i] <= nlstate){ /* a live state */
if(s[m][i]>=1 || s[m][i]==-2 || s[m][i]==-4 || s[m][i]==-5)
mw[++mi][i]=m;
if(m >=lastpass)
else
m++;
}/* end while */
- if (s[m][i] > nlstate){
+ if (s[m][i] > nlstate){ /* In a death state */
mi++; /* Death is another wave */
/* if(mi==0) never been interviewed correctly before death */
/* Only death is a correct wave */
mw[mi][i]=m;
+ }else if (andc[i] != 9999) { /* A death occured after lastpass */
+ m++;
+ mi++;
+ s[m][i]=nlstate+1; /* We are setting the status to the last of non live state */
+ mw[mi][i]=m;
+ nbwarn++;
+ if(firstwo==0){
+ printf("Warning! Death for individual %ld line=%d occurred after last wave %d. Since 0.98r4 we considered a status %d at wave %d\nOthers in log file only\n",num[i],i,lastpass,nlstate+1, m);
+ fprintf(ficlog,"Warning! Death for individual %ld line=%d occurred after last wave %d. Since 0.98r4 we considered a status %d at wave %d\n",num[i],i,lastpass,nlstate+1, m);
+ firstwo=1;
+ }
+ if(firstwo==1){
+ fprintf(ficlog,"Warning! Death for individual %ld line=%d occurred after last wave %d. Since 0.98r4 we considered a status %d at wave %d\n",num[i],i,lastpass,nlstate+1, m);
+ }
}
-
wav[i]=mi;
if(mi==0){
nbwarn++;
}
} /* end mi==0 */
} /* End individuals */
+ /* wav and mw are no more changed */
+
for(i=1; i<=imx; i++){
for(mi=1; mi<wav[i];mi++){
if (stepm <=0)
fprintf(fichtm,"<ul><li><a href='#firstorder'>Result files (first order: no variance)</a>\n \
<li><a href='#secondorder'>Result files (second order (variance)</a>\n \
</ul>");
- fprintf(fichtm,"<ul><li><h4><a name='firstorder'>Result files (first order: no variance)</a></h4>\n \
- - Observed prevalence in each state (during the period defined between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf): <a href=\"%s\">%s</a> (html file) ",
+ fprintf(fichtm,"<ul><li><h4><a name='firstorder'>Result files (first order: no variance)</a></h4>\n");
+ fprintf(fichtm,"<li>- Observed frequency between two states (during the period defined between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf): <a href=\"%s\">%s</a> (html file)<br/>\n",
+ jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,subdirfext3(optionfilefiname,"PHTMFR_",".htm"),subdirfext3(optionfilefiname,"PHTMFR_",".htm"));
+ fprintf(fichtm,"<li> - Observed prevalence in each state (during the period defined between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf): <a href=\"%s\">%s</a> (html file) ",
jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,subdirfext3(optionfilefiname,"PHTM_",".htm"),subdirfext3(optionfilefiname,"PHTM_",".htm"));
fprintf(fichtm,", <a href=\"%s\">%s</a> (text file) <br>\n",subdirf2(fileresu,"P_"),subdirf2(fileresu,"P_"));
fprintf(fichtm,"\
if(prevfcast==1){
/* Projection of prevalence up to period (stable) prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"<br>\n- Projection of cross-sectional prevalence (estimated with cases observed from %
f to %f) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> \
+ fprintf(fichtm,"<br>\n- Projection of cross-sectional prevalence (estimated with cases observed from %
.1f to %.1f) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> \
<img src=\"%s_%d-%d.svg\">", dateprev1, dateprev2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1);
}
}
/* replot exp(p1+p2*x)/(1+exp(p1+p2*x)+exp(p3+p4*x)+exp(p5+p6*x)) t "p12(x)" */
/* fprintf(ficgp,"\nset out \"%s.svg\";",subdirf2(optionfilefiname,"ILK_")); */
fprintf(ficgp,"\nset out \"%s-dest.png\";",subdirf2(optionfilefiname,"ILK_"));
- fprintf(ficgp,"\nset log y;plot \"%s\" u 2:(-$12):5 t \"All sample, transitions colored by destination\" with dots lc variable; set out;\n",subdirf(fileresilk));
+ fprintf(ficgp,"\nset log y;plot \"%s\" u 2:(-$13):6 t \"All sample, transitions colored by destination\" with dots lc variable; set out;\n",subdirf(fileresilk));
fprintf(ficgp,"\nset out \"%s-ori.png\";",subdirf2(optionfilefiname,"ILK_"));
- fprintf(ficgp,"\nset log y;plot \"%s\" u 2:(-$12):4 t \"All sample, transitions colored by origin\" with dots lc variable; set out;\n\n",subdirf(fileresilk));
+ fprintf(ficgp,"\nset log y;plot \"%s\" u 2:(-$13):5 t \"All sample, transitions colored by origin\" with dots lc variable; set out;\n\n",subdirf(fileresilk));
for (i=1; i<= nlstate ; i ++) {
fprintf(ficgp,"\nset out \"%s-p%dj.png\";set ylabel \"Probability for each individual/wave\";",subdirf2(optionfilefiname,"ILK_"),i);
fprintf(ficgp,"unset log;\n# plot weighted, mean weight should have point size of 0.5\n plot \"%s\"",subdirf(fileresilk));
- fprintf(ficgp," u 2:($4 == %d && $5==%d ? $9 : 1/0):($11/4.):5 t \"p%d%d\" with points pointtype 7 ps variable lc variable \\\n",i,1,i,1);
+ fprintf(ficgp," u 2:($5 == %d && $6==%d ? $10 : 1/0):($12/4.):6 t \"p%d%d\" with points pointtype 7 ps variable lc variable \\\n",i,1,i,1);
for (j=2; j<= nlstate+ndeath ; j ++) {
- fprintf(ficgp,",\\\n \"\" u 2:($4 == %d && $5==%d ? $9 : 1/0):($11/4.):5 t \"p%d%d\" with points pointtype 7 ps variable lc variable ",i,j,i,j);
+ fprintf(ficgp,",\\\n \"\" u 2:($5 == %d && $6==%d ? $10 : 1/0):($12/4.):6 t \"p%d%d\" with points pointtype 7 ps variable lc variable ",i,j,i,j);
}
fprintf(ficgp,";\nset out; unset ylabel;\n");
}
in each health status at the date of interview (if between dateprev1 and dateprev2).
We still use firstpass and lastpass as another selection.
*/
+ /* freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint,strstart,\ */
+ /* firstpass, lastpass, stepm, weightopt, model); */
prevalence(probs, ageminpar, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
strcpy(fileresf,"F_");
for (i=1; i<=imx; i++) {
agedc[i]=(moisdc[i]/12.+andc[i])-(moisnais[i]/12.+annais[i]);
for(m=firstpass; (m<= lastpass); m++){
- if(s[m][i] >0 || s[m][i]==-2 || s[m][i]==-4 || s[m][i]==-5){
+ if(s[m][i] >0 || s[m][i]==-1 || s[m][i]==-2 || s[m][i]==-4 || s[m][i]==-5){ /* What if s[m][i]=-1 */
if (s[m][i] >= nlstate+1) {
if(agedc[i]>0){
if((int)moisdc[i]!=99 && (int)andc[i]!=9999){
agev[m][i]=agedc[i];
- /*if(moisdc[i]==99 && andc[i]==9999) s[m][i]=-1;*/
+ /*if(moisdc[i]==99 && andc[i]==9999) s[m][i]=-1;*/
}else {
if ((int)andc[i]!=9999){
nbwarn++;
}
}
} /* agedc > 0 */
- }
+ } /* end if */
else if(s[m][i] !=9){ /* Standard case, age in fractional
years but with the precision of a month */
agev[m][i]=(mint[m][i]/12.+1./24.+anint[m][i])-(moisnais[i]/12.+1./24.+annais[i]);
}
/*agev[m][i]=anint[m][i]-annais[i];*/
/* agev[m][i] = age[i]+2*m;*/
- }
+ } /* en if 9*/
else { /* =9 */
+ /* printf("Debug num[%d]=%ld s[%d][%d]=%d\n",i,num[i], m,i, s[m][i]); */
agev[m][i]=1;
s[m][i]=-1;
}
}
- else /*= 0 Unknown */
+ else if(s[m][i]==0) /*= 0 Unknown */
agev[m][i]=1;
- }
-
+ else{
+ printf("Warning, num[%d]=%ld, s[%d][%d]=%d\n", i, num[i], m, i,s[m][i]);
+ fprintf(ficlog, "Warning, num[%d]=%ld, s[%d][%d]=%d\n", i, num[i], m, i,s[m][i]);
+ agev[m][i]=0;
+ }
+ } /* End for lastpass */
}
+
for (i=1; i<=imx; i++) {
for(m=firstpass; (m<=lastpass); m++){
if (s[m][i] > (nlstate+ndeath)) {
/* */
wav=ivector(1,imx);
- dh=imatrix(1,lastpass-firstpass+1,1,imx);
- bh=imatrix(1,lastpass-firstpass+1,1,imx);
- mw=imatrix(1,lastpass-firstpass+1,1,imx);
+ /* dh=imatrix(1,lastpass-firstpass+1,1,imx); */
+ /* bh=imatrix(1,lastpass-firstpass+1,1,imx); */
+ /* mw=imatrix(1,lastpass-firstpass+1,1,imx); */
+ dh=imatrix(1,lastpass-firstpass+2,1,imx); /* We are adding a wave if status is unknown at last wave but death occurs after last wave.*/
+ bh=imatrix(1,lastpass-firstpass+2,1,imx);
+ mw=imatrix(1,lastpass-firstpass+2,1,imx);
/* Concatenates waves */
+ /* Concatenates waves: wav[i] is the number of effective (useful waves) of individual i.
+ Death is a valid wave (if date is known).
+ mw[mi][i] is the number of (mi=1 to wav[i]) effective wave out of mi of individual i
+ dh[m][i] or dh[mw[mi][i]][i] is the delay between two effective waves m=mw[mi][i]
+ and mw[mi+1][i]. dh depends on stepm.
+ */
+
concatwav(wav, dh, bh, mw, s, agedc, agev, firstpass, lastpass, imx, nlstate, stepm);
/* */
/* Calculates basic frequencies. Computes observed prevalence at single age
and prints on file fileres'p'. */
- freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint,strstart);
+ freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint,strstart,\
+ firstpass, lastpass, stepm, weightopt, model);
fprintf(fichtm,"\n");
fprintf(fichtm,"<br>Total number of observations=%d <br>\n\
/* chdir(path); */
free_ivector(wav,1,imx);
- free_imatrix(dh,1,lastpass-firstpass+1,1,imx);
- free_imatrix(bh,1,lastpass-firstpass+1,1,imx);
- free_imatrix(mw,1,lastpass-firstpass+1,1,imx);
+ free_imatrix(dh,1,lastpass-firstpass+2,1,imx);
+ free_imatrix(bh,1,lastpass-firstpass+2,1,imx);
+ free_imatrix(mw,1,lastpass-firstpass+2,1,imx);
free_lvector(num,1,n);
free_vector(agedc,1,n);
/*free_matrix(covar,0,NCOVMAX,1,n);*/
git://henry.ined.fr / .git / commitdiff