--- imach/src/imach.c 2016/09/15 15:01:13 1.251
+++ imach/src/imach.c 2017/04/04 13:01:16 1.259
@@ -1,6 +1,32 @@
-/* $Id: imach.c,v 1.251 2016/09/15 15:01:13 brouard Exp $
+/* $Id: imach.c,v 1.259 2017/04/04 13:01:16 brouard Exp $
$State: Exp $
$Log: imach.c,v $
+ Revision 1.259 2017/04/04 13:01:16 brouard
+ Summary: Some errors to warnings only if date of death is unknown but status is death we could set to pi3
+
+ Revision 1.258 2017/04/03 10:17:47 brouard
+ Summary: Version 0.99r12
+
+ Some cleanings, conformed with updated documentation.
+
+ Revision 1.257 2017/03/29 16:53:30 brouard
+ Summary: Temp
+
+ Revision 1.256 2017/03/27 05:50:23 brouard
+ Summary: Temporary
+
+ Revision 1.255 2017/03/08 16:02:28 brouard
+ Summary: IMaCh version 0.99r10 bugs in gnuplot fixed
+
+ Revision 1.254 2017/03/08 07:13:00 brouard
+ Summary: Fixing data parameter line
+
+ Revision 1.253 2016/12/15 11:59:41 brouard
+ Summary: 0.99 in progress
+
+ Revision 1.252 2016/09/15 21:15:37 brouard
+ *** empty log message ***
+
Revision 1.251 2016/09/15 15:01:13 brouard
Summary: not working
@@ -128,9 +154,7 @@
Author: Nicolas Brouard
Revision 1.210 2015/11/18 17:41:20 brouard
- Summary: Start working on projected prevalences
-
- Revision 1.209 2015/11/17 22:12:03 brouard
+ Summary: Start working on projected prevalences Revision 1.209 2015/11/17 22:12:03 brouard
Summary: Adding ftolpl parameter
Author: N Brouard
@@ -955,12 +979,12 @@ typedef struct {
#define ODIRSEPARATOR '\\'
#endif
-/* $Id: imach.c,v 1.251 2016/09/15 15:01:13 brouard Exp $ */
+/* $Id: imach.c,v 1.259 2017/04/04 13:01:16 brouard Exp $ */
/* $State: Exp $ */
#include "version.h"
char version[]=__IMACH_VERSION__;
char copyright[]="February 2016,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2018";
-char fullversion[]="$Revision: 1.251 $ $Date: 2016/09/15 15:01:13 $";
+char fullversion[]="$Revision: 1.259 $ $Date: 2017/04/04 13:01:16 $";
char strstart[80];
char optionfilext[10], optionfilefiname[FILENAMELENGTH];
int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */
@@ -1161,6 +1185,7 @@ int *TvarsQind;
#define MAXRESULTLINES 10
int nresult=0;
+int parameterline=0; /* # of the parameter (type) line */
int TKresult[MAXRESULTLINES];
int Tresult[MAXRESULTLINES][NCOVMAX];/* For dummy variable , value (output) */
int Tinvresult[MAXRESULTLINES][NCOVMAX];/* For dummy variable , value (output) */
@@ -2393,8 +2418,8 @@ void powell(double p[], double **xi, int
flatd++;
}
if(flatd >0){
- printf("%d flat directions\n",flatd);
- fprintf(ficlog,"%d flat directions\n",flatd);
+ printf("%d flat directions: ",flatd);
+ fprintf(ficlog,"%d flat directions :",flatd);
for (j=1;j<=n;j++) {
if(flatdir[j]>0){
printf("%d ",j);
@@ -4122,7 +4147,16 @@ void ludcmp(double **a, int n, int *indx
big=0.0;
for (j=1;j<=n;j++)
if ((temp=fabs(a[i][j])) > big) big=temp;
- if (big == 0.0) nrerror("Singular matrix in routine ludcmp");
+ if (big == 0.0){
+ printf(" Singular Hessian matrix at row %d:\n",i);
+ for (j=1;j<=n;j++) {
+ printf(" a[%d][%d]=%f,",i,j,a[i][j]);
+ fprintf(ficlog," a[%d][%d]=%f,",i,j,a[i][j]);
+ }
+ fflush(ficlog);
+ fclose(ficlog);
+ nrerror("Singular matrix in routine ludcmp");
+ }
vv[i]=1.0/big;
}
for (j=1;j<=n;j++) {
@@ -4188,17 +4222,84 @@ void pstamp(FILE *fichier)
fprintf(fichier,"# %s.%s\n#IMaCh version %s, %s\n#%s\n# %s", optionfilefiname,optionfilext,version,copyright, fullversion, strstart);
}
+int linreg(int ifi, int ila, int *no, const double x[], const double y[], double* a, double* b, double* r, double* sa, double * sb) {
+
+ /* y=a+bx regression */
+ double sumx = 0.0; /* sum of x */
+ double sumx2 = 0.0; /* sum of x**2 */
+ double sumxy = 0.0; /* sum of x * y */
+ double sumy = 0.0; /* sum of y */
+ double sumy2 = 0.0; /* sum of y**2 */
+ double sume2; /* sum of square or residuals */
+ double yhat;
+
+ double denom=0;
+ int i;
+ int ne=*no;
+
+ for ( i=ifi, ne=0;i<=ila;i++) {
+ if(!isfinite(x[i]) || !isfinite(y[i])){
+ /* printf(" x[%d]=%f, y[%d]=%f\n",i,x[i],i,y[i]); */
+ continue;
+ }
+ ne=ne+1;
+ sumx += x[i];
+ sumx2 += x[i]*x[i];
+ sumxy += x[i] * y[i];
+ sumy += y[i];
+ sumy2 += y[i]*y[i];
+ denom = (ne * sumx2 - sumx*sumx);
+ /* printf("ne=%d, i=%d,x[%d]=%f, y[%d]=%f sumx=%f, sumx2=%f, sumxy=%f, sumy=%f, sumy2=%f, denom=%f\n",ne,i,i,x[i],i,y[i], sumx, sumx2,sumxy, sumy, sumy2,denom); */
+ }
+
+ denom = (ne * sumx2 - sumx*sumx);
+ if (denom == 0) {
+ // vertical, slope m is infinity
+ *b = INFINITY;
+ *a = 0;
+ if (r) *r = 0;
+ return 1;
+ }
+
+ *b = (ne * sumxy - sumx * sumy) / denom;
+ *a = (sumy * sumx2 - sumx * sumxy) / denom;
+ if (r!=NULL) {
+ *r = (sumxy - sumx * sumy / ne) / /* compute correlation coeff */
+ sqrt((sumx2 - sumx*sumx/ne) *
+ (sumy2 - sumy*sumy/ne));
+ }
+ *no=ne;
+ for ( i=ifi, ne=0;i<=ila;i++) {
+ if(!isfinite(x[i]) || !isfinite(y[i])){
+ /* printf(" x[%d]=%f, y[%d]=%f\n",i,x[i],i,y[i]); */
+ continue;
+ }
+ ne=ne+1;
+ yhat = y[i] - *a -*b* x[i];
+ sume2 += yhat * yhat ;
+
+ denom = (ne * sumx2 - sumx*sumx);
+ /* printf("ne=%d, i=%d,x[%d]=%f, y[%d]=%f sumx=%f, sumx2=%f, sumxy=%f, sumy=%f, sumy2=%f, denom=%f\n",ne,i,i,x[i],i,y[i], sumx, sumx2,sumxy, sumy, sumy2,denom); */
+ }
+ *sb = sqrt(sume2/(ne-2)/(sumx2 - sumx * sumx /ne));
+ *sa= *sb * sqrt(sumx2/ne);
+
+ return 0;
+}
+
/************ Frequencies ********************/
void freqsummary(char fileres[], double p[], double pstart[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, \
int *Tvaraff, int *invalidvarcomb, int **nbcode, int *ncodemax,double **mint,double **anint, char strstart[], \
int firstpass, int lastpass, int stepm, int weightopt, char model[])
{ /* Some frequencies as well as proposing some starting values */
- int i, m, jk, j1, bool, z1,j, k, iv, jj=0;
+ int i, m, jk, j1, bool, z1,j, nj, nl, k, iv, jj=0;
int iind=0, iage=0;
int mi; /* Effective wave */
int first;
double ***freq; /* Frequencies */
+ double *x, *y, a,b,r, sa, sb; /* for regression, y=b+m*x and r is the correlation coefficient */
+ int no;
double *meanq;
double **meanqt;
double *pp, **prop, *posprop, *pospropt;
@@ -4251,6 +4352,8 @@ Title=%s
Datafile=%s Firstpass=%d La
}
fprintf(ficresphtmfr,"Current page is file %s
\n\nFrequencies of all effective transitions of the model, by age at begin of transition, and covariate value at the begin of transition (if the covariate is a varying covariate)
Unknown status is -1
\n",fileresphtmfr, fileresphtmfr);
+ y= vector(iagemin-AGEMARGE,iagemax+4+AGEMARGE);
+ x= vector(iagemin-AGEMARGE,iagemax+4+AGEMARGE);
freq= ma3x(-5,nlstate+ndeath,-5,nlstate+ndeath,iagemin-AGEMARGE,iagemax+4+AGEMARGE);
j1=0;
@@ -4268,7 +4371,15 @@ Title=%s
Datafile=%s Firstpass=%d La
dateintsum=0;
k2cpt=0;
- for (j = 0; j <= cptcoveff; j+=cptcoveff){ /* j= 0 constant model */
+ if(cptcoveff == 0 )
+ nl=1; /* Constant model only */
+ else
+ nl=2;
+ for (nj = 1; nj <= nl; nj++){ /* nj= 1 constant model, nl number of loops. */
+ if(nj==1)
+ j=0; /* First pass for the constant */
+ else
+ j=cptcoveff; /* Other passes for the covariate values */
first=1;
for (j1 = 1; j1 <= (int) pow(2,j); j1++){ /* Loop on covariates combination in order of model, excluding quantitatives, V4=0, V3=0 for example, fixed or varying covariates */
posproptt=0.;
@@ -4514,6 +4625,7 @@ Title=%s
Datafile=%s Firstpass=%d La
if(first==1){
printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]);
}
+ /* printf(" %d%d=%.0f",jk,m,freq[jk][m][iage]); */
fprintf(ficlog," %d%d=%.0f",jk,m,freq[jk][m][iage]);
}
if(jk!=0 && m!=0)
@@ -4560,22 +4672,33 @@ Title=%s
Datafile=%s Firstpass=%d La
for(i=1,jk=1; i <=nlstate; i++){
for(k=1; k <=(nlstate+ndeath); k++){
if (k != i) {
- printf("%d%d ",i,k);
- fprintf(ficlog,"%d%d ",i,k);
for(jj=1; jj <=ncovmodel; jj++){ /* For counting jk */
- if(jj==1){ /* Constant case */
+ if(jj==1){ /* Constant case (in fact cste + age) */
if(j1==1){ /* All dummy covariates to zero */
freq[i][k][iagemax+4]=freq[i][k][iagemax+3]; /* Stores case 0 0 0 */
freq[i][i][iagemax+4]=freq[i][i][iagemax+3]; /* Stores case 0 0 0 */
+ printf("%d%d ",i,k);
+ fprintf(ficlog,"%d%d ",i,k);
+ printf("%12.7f ln(%.0f/%.0f)= %f, OR=%f sd=%f \n",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]),freq[i][k][iagemax+3]/freq[i][i][iagemax+3], sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3]));
+ fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f \n",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]));
+ pstart[jk]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]);
}
- printf("%12.7f ln(%.0f/%.0f)= %f, OR=%f sd=%f \n",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]),freq[i][k][iagemax+3]/freq[i][i][iagemax+3], sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3]));
- fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f \n",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]));
- pstart[jk]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]);
- }else if( (log(j1-1)/log(2)+1 == jj -2 -nagesqr) && Dummy[jj-2-nagesqr]==0){ /* We want only if the position, jj, in model corresponds to unique covariate equal to 1 in j1 combination */
+ }else if((j1==1) && (jj==2 || nagesqr==1)){ /* age or age*age parameter without covariate V4*age (to be done later) */
+ for(iage=iagemin; iage <= iagemax+3; iage++){
+ x[iage]= (double)iage;
+ y[iage]= log(freq[i][k][iage]/freq[i][i][iage]);
+ /* printf("i=%d, k=%d, jk=%d, j1=%d, jj=%d, y[%d]=%f\n",i,k,jk,j1,jj, iage, y[iage]); */
+ }
+ linreg(iagemin,iagemax,&no,x,y,&a,&b,&r, &sa, &sb ); /* y= a+b*x with standard errors */
+ pstart[jk]=b;
+ pstart[jk-1]=a;
+ }else if( j1!=1 && (j1==2 || (log(j1-1.)/log(2.)-(int)(log(j1-1.)/log(2.))) <0.010) && ( TvarsDind[(int)(log(j1-1.)/log(2.))+1]+2+nagesqr == jj) && Dummy[jj-2-nagesqr]==0){ /* We want only if the position, jj, in model corresponds to unique covariate equal to 1 in j1 combination */
+ printf("j1=%d, jj=%d, (int)(log(j1-1.)/log(2.))+1=%d, TvarsDind[(int)(log(j1-1.)/log(2.))+1]=%d\n",j1, jj,(int)(log(j1-1.)/log(2.))+1,TvarsDind[(int)(log(j1-1.)/log(2.))+1]);
+ printf("j1=%d, jj=%d, (log(j1-1.)/log(2.))+1=%f, TvarsDind[(int)(log(j1-1.)/log(2.))+1]=%d\n",j1, jj,(log(j1-1.)/log(2.))+1,TvarsDind[(int)(log(j1-1.)/log(2.))+1]);
pstart[jk]= log((freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4]));
+ printf("%d%d ",i,k);
+ fprintf(ficlog,"%d%d ",i,k);
printf("jk=%d,i=%d,k=%d,p[%d]=%12.7f ln((%.0f/%.0f)/(%.0f/%.0f))= %f, OR=%f sd=%f \n",jk,i,k,jk,p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3],freq[i][k][iagemax+4],freq[i][i][iagemax+4], log((freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4])),(freq[i][k][iagemax+3]/freq[i][i][iagemax+3])/(freq[i][k][iagemax+4]/freq[i][i][iagemax+4]), sqrt(1/freq[i][k][iagemax+3]+1/freq[i][i][iagemax+3]+1/freq[i][k][iagemax+4]+1/freq[i][i][iagemax+4]));
- }else if(jj==2 || nagesqr==1){ /* age or age*age parameter */
- ;
}else{ /* Other cases, like quantitative fixed or varying covariates */
;
}
@@ -4583,8 +4706,6 @@ Title=%s
Datafile=%s Firstpass=%d La
/* fprintf(ficlog,"%12.7f )", param[i][jj][k]); */
jk++;
} /* end jj */
- printf("\n");
- fprintf(ficlog,"\n");
} /* end k!= i */
} /* end k */
} /* end i, jk */
@@ -4599,7 +4720,9 @@ Title=%s
Datafile=%s Firstpass=%d La
printf("%d%d ",i,k);
fprintf(ficlog,"%d%d ",i,k);
for(jj=1; jj <=ncovmodel; jj++){
- if(jj==1){
+ pstart[jk]=p[jk]; /* Setting pstart to p values by default */
+ if(jj==1){ /* Age has to be done */
+ pstart[jk]= log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]);
printf("%12.7f ln(%.0f/%.0f)= %12.7f ",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]));
fprintf(ficlog,"%12.7f ln(%.0f/%.0f)= %12.7f ",p[jk],freq[i][k][iagemax+3],freq[i][i][iagemax+3], log(freq[i][k][iagemax+3]/freq[i][i][iagemax+3]));
}
@@ -4630,6 +4753,28 @@ Title=%s
Datafile=%s Firstpass=%d La
fprintf(ficlog,"\n");
} /* end j=0 */
} /* end j */
+
+ if(mle == -2){ /* We want to use these values as starting values */
+ for(i=1, jk=1; i <=nlstate; i++){
+ for(j=1; j <=nlstate+ndeath; j++){
+ if(j!=i){
+ /*ca[0]= k+'a'-1;ca[1]='\0';*/
+ printf("%1d%1d",i,j);
+ fprintf(ficparo,"%1d%1d",i,j);
+ for(k=1; k<=ncovmodel;k++){
+ /* printf(" %lf",param[i][j][k]); */
+ /* fprintf(ficparo," %lf",param[i][j][k]); */
+ p[jk]=pstart[jk];
+ printf(" %f ",pstart[jk]);
+ fprintf(ficparo," %f ",pstart[jk]);
+ jk++;
+ }
+ printf("\n");
+ fprintf(ficparo,"\n");
+ }
+ }
+ }
+ } /* end mle=-2 */
dateintmean=dateintsum/k2cpt;
fclose(ficresp);
@@ -4637,6 +4782,8 @@ Title=%s
Datafile=%s Firstpass=%d La
fclose(ficresphtmfr);
free_vector(meanq,1,nqfveff);
free_matrix(meanqt,1,lastpass,1,nqtveff);
+ free_vector(x, iagemin-AGEMARGE, iagemax+4+AGEMARGE);
+ free_vector(y, iagemin-AGEMARGE, iagemax+4+AGEMARGE);
free_ma3x(freq,-5,nlstate+ndeath,-5,nlstate+ndeath, iagemin-AGEMARGE, iagemax+4+AGEMARGE);
free_vector(pospropt,1,nlstate);
free_vector(posprop,1,nlstate);
@@ -4824,13 +4971,9 @@ void concatwav(int wav[], int **dh, int
/* if(mi==0) never been interviewed correctly before death */
/* Only death is a correct wave */
mw[mi][i]=m;
- }
+ } /* else not in a death state */
#ifndef DISPATCHINGKNOWNDEATHAFTERLASTWAVE
- else if ((int) andc[i] != 9999) { /* Status is negative. A death occured after lastpass, we can't take it into account because of potential bias */
- /* m++; */
- /* mi++; */
- /* s[m][i]=nlstate+1; /\* We are setting the status to the last of non live state *\/ */
- /* mw[mi][i]=m; */
+ else if ((int) andc[i] != 9999) { /* Date of death is known */
if ((int)anint[m][i]!= 9999) { /* date of last interview is known */
if((andc[i]+moisdc[i]/12.) <=(anint[m][i]+mint[m][i]/12.)){ /* death occured before last wave and status should have been death instead of -1 */
nbwarn++;
@@ -4843,12 +4986,12 @@ void concatwav(int wav[], int **dh, int
}else{ /* Death occured afer last wave potential bias */
nberr++;
if(firstwo==0){
- printf("Error! Death for individual %ld line=%d occurred at %d/%d after last wave %d interviewed at %d/%d. Potential bias if other individuals are still alive at this date but ignored. This case (%d)/wave (%d) is skipped, no contribution to likelihood.\nOthers in log file only\n",num[i],i,(int) moisdc[i], (int) andc[i], lastpass,(int)mint[m][i],(int)anint[m][i], i,m );
+ printf("Error! Death for individual %ld line=%d occurred at %d/%d after last wave %d interviewed at %d/%d. Potential bias if other individuals are still alive at this date but ignored. This case (%d)/wave (%d) is skipped, no contribution to likelihood. Please add a new fictive wave at the date of last vital status scan, with a dead status or alive but unknown state status (-1). See documentation\nOthers in log file only\n",num[i],i,(int) moisdc[i], (int) andc[i], lastpass,(int)mint[m][i],(int)anint[m][i], i,m );
firstwo=1;
}
- fprintf(ficlog,"Error! Death for individual %ld line=%d occurred at %d/%d after last wave %d interviewed at %d/%d. Potential bias if other individuals are still alive at this date but ignored. This case (%d)/wave (%d) is skipped, no contribution to likelihood.\n",num[i],i,(int) moisdc[i], (int) andc[i], lastpass,(int)mint[m][i],(int)anint[m][i], i,m );
+ fprintf(ficlog,"Error! Death for individual %ld line=%d occurred at %d/%d after last wave %d interviewed at %d/%d. Potential bias if other individuals are still alive at this date but ignored. This case (%d)/wave (%d) is skipped, no contribution to likelihood. Please add a new fictive wave at the date of last vital status scan, with a dead status or alive but unknown state status (-1). See documentation\n\n",num[i],i,(int) moisdc[i], (int) andc[i], lastpass,(int)mint[m][i],(int)anint[m][i], i,m );
}
- }else{ /* end date of interview is known */
+ }else{ /* if date of interview is unknown */
/* death is known but not confirmed by death status at any wave */
if(firstfour==0){
printf("Error! Death for individual %ld line=%d occurred %d/%d but not confirmed by any death status for any wave, including last wave %d at unknown date %d/%d. Potential bias if other individuals are still alive at this date but ignored. This case (%d)/wave (%d) is skipped, no contribution to likelihood.\nOthers in log file only\n",num[i],i,(int) moisdc[i], (int) andc[i], lastpass,(int)mint[m][i],(int)anint[m][i], i,m );
@@ -6241,7 +6384,7 @@ To be simple, these graphs help to under
void printinghtml(char fileresu[], char title[], char datafile[], int firstpass, \
int lastpass, int stepm, int weightopt, char model[],\
int imx,int jmin, int jmax, double jmeanint,char rfileres[],\
- int popforecast, int prevfcast, int backcast, int estepm , \
+ int popforecast, int mobilav, int prevfcast, int mobilavproj, int backcast, int estepm , \
double jprev1, double mprev1,double anprev1, double dateprev1, \
double jprev2, double mprev2,double anprev2, double dateprev2){
int jj1, k1, i1, cpt, k4, nres;
@@ -6288,7 +6431,7 @@ void printinghtml(char fileresu[], char
for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k1=1; k1<=m;k1++){ /* For each combination of covariate */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
/* for(i1=1; i1<=ncodemax[k1];i1++){ */
@@ -6315,7 +6458,7 @@ void printinghtml(char fileresu[], char
}
}
/* aij, bij */
- fprintf(fichtm,"
- Logit model (yours is: 1+age+%s), for example: logit(pij)=log(pij/pii)= aij+ bij age + V1 age + etc. as a function of age: %s_%d-1-%d.svg
\
+ fprintf(fichtm,"
- Logit model (yours is: logit(pij)=log(pij/pii)= aij+ bij age+%s) as a function of age: %s_%d-1-%d.svg
\
![](\"%s_%d-1-%d.svg\")
",model,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres);
/* Pij */
fprintf(fichtm,"
\n- Pij or conditional probabilities to be observed in state j being in state i, %d (stepm) months before: %s_%d-2-%d.svg
\
@@ -6339,21 +6482,21 @@ divided by h: hPij
}
/* Period (stable) prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\n- Convergence to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+ fprintf(fichtm,"
\n- Convergence to period (stable) prevalence in state %d. Or probability to be in state %d some years earlier, knowing that we will be in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
![](\"%s_%d-%d-%d.svg\")
", cpt, cpt, nlstate, subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres);
}
if(backcast==1){
/* Period (stable) back prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\n- Convergence to period (stable) back prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+ fprintf(fichtm,"
\n- Convergence to mixed (stable) back prevalence in state %d. Or probability to be in state %d at a younger age, knowing that we will be in state (1 to %d) at different older ages. %s_%d-%d-%d.svg
\
", cpt, cpt, nlstate, subdirf2(optionfilefiname,"PB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PB_"),cpt,k1,nres,subdirf2(optionfilefiname,"PB_"),cpt,k1,nres);
}
}
if(prevfcast==1){
/* Projection of prevalence up to period (stable) prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
- fprintf(fichtm,"
\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. %s_%d-%d-%d.svg
\
-", dateprev1, dateprev2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres);
+ fprintf(fichtm,"
\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f and mobil_average=%d) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. %s_%d-%d-%d.svg
\
+", dateprev1, dateprev2, mobilavproj, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres);
}
}
@@ -6420,7 +6563,7 @@ See page 'Matrix of variance-covariance
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
for(k1=1; k1<=m;k1++){
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
/* for(i1=1; i1<=ncodemax[k1];i1++){ */
jj1++;
@@ -6441,9 +6584,9 @@ See page 'Matrix of variance-covariance
}
}
for(cpt=1; cpt<=nlstate;cpt++) {
- fprintf(fichtm,"\n
- Observed (cross-sectional) and period (incidence based) \
+ fprintf(fichtm,"\n
- Observed (cross-sectional with mov_average=%d) and period (incidence based) \
prevalence (with 95%% confidence interval) in state (%d): %s_%d-%d-%d.svg\n
\
-",cpt,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres);
+",mobilav,cpt,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres,subdirf2(optionfilefiname,"V_"),cpt,k1,nres);
}
fprintf(fichtm,"\n
- Total life expectancy by age and \
health expectancies in states (1) and (2). If popbased=1 the smooth (due to the model) \
@@ -6516,7 +6659,7 @@ void printinggnuplot(char fileresu[], ch
for (k1=1; k1<= m ; k1 ++){ /* For each valid combination of covariate */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
/* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
/* We are interested in selected combination by the resultline */
/* printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt); */
@@ -6598,7 +6741,7 @@ void printinggnuplot(char fileresu[], ch
/*2 eme*/
for (k1=1; k1<= m ; k1 ++){
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files ");
for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
@@ -6658,7 +6801,7 @@ void printinggnuplot(char fileresu[], ch
/*3eme*/
for (k1=1; k1<= m ; k1 ++){
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
for (cpt=1; cpt<= nlstate ; cpt ++) {
@@ -6707,7 +6850,7 @@ plot [%.f:%.f] \"%s\" every :::%d::%d u
/* Survival functions (period) from state i in state j by initial state i */
for (k1=1; k1<=m; k1++){ /* For each covariate and each value */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state cpt*/
fprintf(ficgp,"\n#\n#\n# Survival functions in state j : 'LIJ_' files, cov=%d state=%d",k1, cpt);
@@ -6753,7 +6896,7 @@ set ter svg size 640, 480\nunset log y\n
/* Survival functions (period) from state i in state j by final state j */
for (k1=1; k1<= m ; k1++){ /* For each covariate combination if any */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state */
fprintf(ficgp,"\n#\n#\n# Survival functions in state j and all livestates from state i by final state j: 'lij' files, cov=%d state=%d",k1, cpt);
@@ -6807,9 +6950,9 @@ set ter svg size 640, 480\nunset log y\n
/* CV preval stable (period) for each covariate */
for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
- for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
+ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state of arrival */
fprintf(ficgp,"\n#\n#\n#CV preval stable (period): 'pij' files, covariatecombination#=%d state=%d",k1, cpt);
for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
@@ -6833,12 +6976,12 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\
set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
k=3; /* Offset */
- for (i=1; i<= nlstate ; i ++){
+ for (i=1; i<= nlstate ; i ++){ /* State of origin */
if(i==1)
fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
else
fprintf(ficgp,", '' ");
- l=(nlstate+ndeath)*(i-1)+1;
+ l=(nlstate+ndeath)*(i-1)+1; /* 1, 1+ nlstate+ndeath, 1+2*(nlstate+ndeath) */
fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l);
for (j=2; j<= nlstate ; j ++)
fprintf(ficgp,"+$%d",k+l+j-1);
@@ -6854,10 +6997,10 @@ set ter svg size 640, 480\nunset log y\n
/* CV back preval stable (period) for each covariate */
for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
- for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
- fprintf(ficgp,"\n#\n#\n#CV Back preval stable (period): 'pij' files, covariatecombination#=%d state=%d",k1, cpt);
+ for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life ending state */
+ fprintf(ficgp,"\n#\n#\n#CV Back preval stable (period): 'pijb' files, covariatecombination#=%d state=%d",k1, cpt);
for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */
lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
/* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */
@@ -6879,16 +7022,16 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\
set ter svg size 640, 480\nunset log y\nplot [%.f:%.f] ", ageminpar, agemaxpar);
k=3; /* Offset */
- for (i=1; i<= nlstate ; i ++){
+ for (i=1; i<= nlstate ; i ++){ /* State of origin */
if(i==1)
fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJB_"));
else
fprintf(ficgp,", '' ");
/* l=(nlstate+ndeath)*(i-1)+1; */
- l=(nlstate+ndeath)*(cpt-1)+1;
+ l=(nlstate+ndeath)*(cpt-1)+1; /* fixed for i; cpt=1 1, cpt=2 1+ nlstate+ndeath, 1+2*(nlstate+ndeath) */
/* fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l); /\* a vérifier *\/ */
/* fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l+(cpt-1)+i-1); /\* a vérifier *\/ */
- fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d",k1,k+l+(cpt-1)+i-1); /* a vérifier */
+ fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d",k1,k+l+i-1); /* To be verified */
/* for (j=2; j<= nlstate ; j ++) */
/* fprintf(ficgp,"+$%d",k+l+j-1); */
/* /\* fprintf(ficgp,"+$%d",k+l+j-1); *\/ */
@@ -6905,7 +7048,7 @@ set ter svg size 640, 480\nunset log y\n
for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k1)
+ if(m != 1 && TKresult[nres]!= k1)
continue;
for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */
fprintf(ficgp,"\n#\n#\n#Projection of prevalence to stable (period): 'PROJ_' files, covariatecombination#=%d state=%d",k1, cpt);
@@ -7036,7 +7179,7 @@ set ter svg size 640, 480\nunset log y\n
fprintf(ficgp,"# jk=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */
for(jk=1; jk <=m; jk++) /* For each combination of covariate */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= jk)
+ if(m != 1 && TKresult[nres]!= jk)
continue;
fprintf(ficgp,"# Combination of dummy jk=%d and ",jk);
for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
@@ -7398,7 +7541,7 @@ set ter svg size 640, 480\nunset log y\n
/* if (h==(int)(YEARM*yearp)){ */
for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k=1; k<=i1;k++){
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
if(invalidvarcomb[k]){
printf("\nCombination (%d) projection ignored because no cases \n",k);
@@ -8349,6 +8492,11 @@ int decoderesult ( char resultline[], in
if (strlen(resultsav) >1){
j=nbocc(resultsav,'='); /**< j=Number of covariate values'=' */
}
+ if(j == 0){ /* Resultline but no = */
+ TKresult[nres]=0; /* Combination for the nresult and the model */
+ return (0);
+ }
+
if( j != cptcovs ){ /* Be careful if a variable is in a product but not single */
printf("ERROR: the number of variable in the resultline, %d, differs from the number of variable used in the model line, %d.\n",j, cptcovs);
fprintf(ficlog,"ERROR: the number of variable in the resultline, %d, differs from the number of variable used in the model line, %d.\n",j, cptcovs);
@@ -8960,19 +9108,19 @@ int calandcheckages(int imx, int maxwav,
s[m][i]=-1;
}
if((int)moisdc[i]==99 && (int)andc[i]==9999 && s[m][i]>nlstate){
- *nberr = *nberr + 1;
+ *nbwarn = *nbwarn + 1;
if(firstone == 0){
firstone=1;
- printf("Error! Date of death (month %2d and year %4d) of individual %ld on line %d was unknown, you must set an arbitrary year of death or he/she is skipped and results can be biased (%d) because status is a death state %d at wave %d. Wave dropped.\nOther similar cases in log file\n",(int)moisdc[i],(int)andc[i],num[i],i, *nberr,s[m][i],m);
+ printf("Warning (#%d)! Date of death (month %2d and year %4d) of individual %ld on line %d was unknown but status is a death state %d at wave %d. If you don't know the vital status, please enter -2. If he/she is still alive but don't know the state, please code with '-1 or '.'. Here, we believe in a death.\nOther similar cases in log file\n", *nbwarn,(int)moisdc[i],(int)andc[i],num[i],i,s[m][i],m);
}
- fprintf(ficlog,"Error! Date of death (month %2d and year %4d) of individual %ld on line %d was unknown, you must set an arbitrary year of death or he/she is skipped and results can be biased (%d) because status is a death state %d at wave %d. Wave dropped.\n",(int)moisdc[i],(int)andc[i],num[i],i, *nberr,s[m][i],m);
- s[m][i]=-1;
+ fprintf(ficlog,"Warning (#%d)! Date of death (month %2d and year %4d) of individual %ld on line %d was unknown but status is a death state %d at wave %d. If you don't know the vital status, please enter -2. If he/she is still alive but don't know the state, please code with '-1 or '.'. Here, we believe in a death.\nOther similar cases in log file\n", *nbwarn,(int)moisdc[i],(int)andc[i],num[i],i,s[m][i],m);
+ /* s[m][i]=-1; */ /* Keeping the death status */
}
if((int)moisdc[i]==99 && (int)andc[i]!=9999 && s[m][i]>nlstate){
(*nberr)++;
- printf("Error! Month of death of individual %ld on line %d was unknown %2d, you should set it otherwise the information on the death is skipped and results are biased.\n",num[i],i,(int)moisdc[i]);
- fprintf(ficlog,"Error! Month of death of individual %ld on line %d was unknown %f, you should set it otherwise the information on the death is skipped and results are biased.\n",num[i],i,moisdc[i]);
- s[m][i]=-1; /* We prefer to skip it (and to skip it in version 0.8a1 too */
+ printf("Error (#%d)! Month of death of individual %ld on line %d was unknown (%2d) (year of death is %4d) and status is a death state %d at wave %d. Please impute an arbitrary (or not) month and rerun. Currently this transition to death will be skipped (status is set to -2).\nOther similar cases in log file\n", *nberr, num[i],i,(int)moisdc[i],(int)andc[i],s[m][i],m);
+ fprintf(ficlog,"Error (#%d)! Month of death of individual %ld on line %d was unknown (%2d) (year of death is %4d) and status is a death state %d at wave %d. Please impute an arbitrary (or not) month and rerun. Currently this transition to death will be skipped (status is set to -2).\nOther similar cases in log file\n", *nberr, num[i],i,(int)moisdc[i],(int)andc[i],s[m][i],m);
+ s[m][i]=-2; /* We prefer to skip it (and to skip it in version 0.8a1 too */
}
}
}
@@ -9285,7 +9433,7 @@ int prevalence_limit(double *p, double *
for(k=1; k<=i1;k++){ /* For each combination k of dummy covariates in the model */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
/* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
@@ -9382,7 +9530,7 @@ int back_prevalence_limit(double *p, dou
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
//printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
fprintf(ficresplb,"#******");
@@ -9441,6 +9589,7 @@ int back_prevalence_limit(double *p, dou
fprintf(ficresplb," %.3f %d\n", tot, *ncvyearp);
} /* Age */
/* was end of cptcod */
+ /*fprintf(ficresplb,"\n");*/ /* Seems to be necessary for gnuplot only if two result lines and no covariate. */
} /* end of any combination */
} /* end of nres */
/* hBijx(p, bage, fage); */
@@ -9485,7 +9634,7 @@ int hPijx(double *p, int bage, int fage)
/* k=k+1; */
for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k=1; k<=i1;k++){
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
fprintf(ficrespij,"\n#****** ");
for(j=1;j<=cptcoveff;j++)
@@ -9557,14 +9706,14 @@ int hPijx(double *p, int bage, int fage)
/* hstepm=1; aff par mois*/
pstamp(ficrespijb);
- fprintf(ficrespijb,"#****** h Pij x Back Probability to be in state i at age x-h being in j at x ");
+ fprintf(ficrespijb,"#****** h Bij x Back probability to be in state i at age x-h being in j at x: B1j+B2j+...=1 ");
i1= pow(2,cptcoveff);
/* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
/* /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */
/* k=k+1; */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
fprintf(ficrespijb,"\n#****** ");
for(j=1;j<=cptcoveff;j++)
@@ -9591,7 +9740,7 @@ int hPijx(double *p, int bage, int fage)
/* hbxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k); */
hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm, k);
/* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */
- fprintf(ficrespijb,"# Cov Agex agex-h hpijx with i,j=");
+ fprintf(ficrespijb,"# Cov Agex agex-h hbijx with i,j=");
for(i=1; i<=nlstate;i++)
for(j=1; j<=nlstate+ndeath;j++)
fprintf(ficrespijb," %1d-%1d",i,j);
@@ -9637,6 +9786,7 @@ int main(int argc, char *argv[])
int NDIM=2;
int vpopbased=0;
int nres=0;
+ int endishere=0;
char ca[32], cb[32];
/* FILE *fichtm; *//* Html File */
@@ -10007,11 +10157,6 @@ int main(int argc, char *argv[])
fclose (ficlog);
goto end;
exit(0);
- } else if(mle==-2) { /* Guessing from means */
- prwizard(ncovmodel, nlstate, ndeath, model, ficparo);
- printf(" You chose mle=-2, look at file %s for a template of covariance matrix \n",filereso);
- fprintf(ficlog," You chose mle=-2, look at file %s for a template of covariance matrix \n",filereso);
-
} else if(mle==-5) { /* Main Wizard */
prwizard(ncovmodel, nlstate, ndeath, model, ficparo);
printf(" You chose mle=-3, look at file %s for a template of covariance matrix \n",filereso);
@@ -10997,70 +11142,31 @@ Please run with mle=-1 to get a correct
fprintf(ficparo,"agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f estepm=%d, ftolpl=%e\n",ageminpar,agemaxpar,bage,fage, estepm, ftolpl);
/* Other stuffs, more or less useful */
- while((c=getc(ficpar))=='#' && c!= EOF){
- ungetc(c,ficpar);
- fgets(line, MAXLINE, ficpar);
- fputs(line,stdout);
- fputs(line,ficparo);
- }
- ungetc(c,ficpar);
-
- fscanf(ficpar,"begin-prev-date=%lf/%lf/%lf end-prev-date=%lf/%lf/%lf mov_average=%d\n",&jprev1, &mprev1,&anprev1,&jprev2, &mprev2,&anprev2,&mobilav);
- fprintf(ficparo,"begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
- fprintf(ficres,"begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
- printf("begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
- fprintf(ficlog,"begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
-
- while((c=getc(ficpar))=='#' && c!= EOF){
- ungetc(c,ficpar);
- fgets(line, MAXLINE, ficpar);
- fputs(line,stdout);
- fputs(line,ficparo);
- }
- ungetc(c,ficpar);
-
-
- dateprev1=anprev1+(mprev1-1)/12.+(jprev1-1)/365.;
- dateprev2=anprev2+(mprev2-1)/12.+(jprev2-1)/365.;
-
- fscanf(ficpar,"pop_based=%d\n",&popbased);
- fprintf(ficlog,"pop_based=%d\n",popbased);
- fprintf(ficparo,"pop_based=%d\n",popbased);
- fprintf(ficres,"pop_based=%d\n",popbased);
-
- while((c=getc(ficpar))=='#' && c!= EOF){
- ungetc(c,ficpar);
- fgets(line, MAXLINE, ficpar);
- fputs(line,stdout);
- fputs(line,ficres);
- fputs(line,ficparo);
+ while(fgets(line, MAXLINE, ficpar)) {
+ /* If line starts with a # it is a comment */
+ if (line[0] == '#') {
+ numlinepar++;
+ fputs(line,stdout);
+ fputs(line,ficparo);
+ fputs(line,ficlog);
+ continue;
+ }else
+ break;
}
- ungetc(c,ficpar);
-
- fscanf(ficpar,"prevforecast=%d starting-proj-date=%lf/%lf/%lf final-proj-date=%lf/%lf/%lf mobil_average=%d\n",&prevfcast,&jproj1,&mproj1,&anproj1,&jproj2,&mproj2,&anproj2,&mobilavproj);
- fprintf(ficparo,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
- printf("prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
- fprintf(ficlog,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
- fprintf(ficres,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
- /* day and month of proj2 are not used but only year anproj2.*/
-
- while((c=getc(ficpar))=='#' && c!= EOF){
- ungetc(c,ficpar);
- fgets(line, MAXLINE, ficpar);
- fputs(line,stdout);
- fputs(line,ficparo);
- fputs(line,ficres);
+
+ if((num_filled=sscanf(line,"begin-prev-date=%lf/%lf/%lf end-prev-date=%lf/%lf/%lf mov_average=%d\n",&jprev1, &mprev1,&anprev1,&jprev2, &mprev2,&anprev2,&mobilav)) !=EOF){
+
+ if (num_filled != 7) {
+ printf("Error: Not 7 (data)parameters in line but %d, for example:begin-prev-date=1/1/1990 end-prev-date=1/6/2004 mov_average=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ fprintf(ficlog,"Error: Not 7 (data)parameters in line but %d, for example:begin-prev-date=1/1/1990 end-prev-date=1/6/2004 mov_average=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ goto end;
+ }
+ printf("begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
+ fprintf(ficparo,"begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
+ fprintf(ficres,"begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
+ fprintf(ficlog,"begin-prev-date=%.lf/%.lf/%.lf end-prev-date=%.lf/%.lf/%.lf mov_average=%d\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,mobilav);
}
- ungetc(c,ficpar);
-
- fscanf(ficpar,"backcast=%d starting-back-date=%lf/%lf/%lf final-back-date=%lf/%lf/%lf mobil_average=%d\n",&backcast,&jback1,&mback1,&anback1,&jback2,&mback2,&anback2,&mobilavproj);
- fprintf(ficparo,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- fprintf(ficlog,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- fprintf(ficres,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
- /* day and month of proj2 are not used but only year anproj2.*/
-
- /* Results */
- nresult=0;
+
while(fgets(line, MAXLINE, ficpar)) {
/* If line starts with a # it is a comment */
if (line[0] == '#') {
@@ -11068,49 +11174,114 @@ Please run with mle=-1 to get a correct
fputs(line,stdout);
fputs(line,ficparo);
fputs(line,ficlog);
- fputs(line,ficres);
continue;
}else
break;
}
- if (!feof(ficpar))
- while((num_filled=sscanf(line,"result:%[^\n]\n",resultline)) !=EOF){
- if (num_filled == 0){
- resultline[0]='\0';
- break;
- } else if (num_filled != 1){
- printf("ERROR %d: result line should be at minimum 'result=' %s\n",num_filled, line);
- }
- nresult++; /* Sum of resultlines */
- printf("Result %d: result=%s\n",nresult, resultline);
- if(nresult > MAXRESULTLINES){
- printf("ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
- fprintf(ficlog,"ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
+
+
+ dateprev1=anprev1+(mprev1-1)/12.+(jprev1-1)/365.;
+ dateprev2=anprev2+(mprev2-1)/12.+(jprev2-1)/365.;
+
+ if((num_filled=sscanf(line,"pop_based=%d\n",&popbased)) !=EOF){
+ if (num_filled != 1) {
+ printf("Error: Not 1 (data)parameters in line but %d, for example:pop_based=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ fprintf(ficlog,"Error: Not 1 (data)parameters in line but %d, for example: pop_based=1\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
goto end;
}
- decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
- fprintf(ficparo,"result: %s\n",resultline);
- fprintf(ficres,"result: %s\n",resultline);
- fprintf(ficlog,"result: %s\n",resultline);
- while(fgets(line, MAXLINE, ficpar)) {
+ printf("pop_based=%d\n",popbased);
+ fprintf(ficlog,"pop_based=%d\n",popbased);
+ fprintf(ficparo,"pop_based=%d\n",popbased);
+ fprintf(ficres,"pop_based=%d\n",popbased);
+ }
+
+ /* Results */
+ nresult=0;
+ do{
+ if(!fgets(line, MAXLINE, ficpar)){
+ endishere=1;
+ parameterline=14;
+ }else if (line[0] == '#') {
/* If line starts with a # it is a comment */
- if (line[0] == '#') {
- numlinepar++;
- fputs(line,stdout);
- fputs(line,ficparo);
- fputs(line,ficres);
- fputs(line,ficlog);
- continue;
- }else
- break;
+ numlinepar++;
+ fputs(line,stdout);
+ fputs(line,ficparo);
+ fputs(line,ficlog);
+ continue;
+ }else if(sscanf(line,"prevforecast=%[^\n]\n",modeltemp))
+ parameterline=11;
+ else if(sscanf(line,"backcast=%[^\n]\n",modeltemp))
+ parameterline=12;
+ else if(sscanf(line,"result:%[^\n]\n",modeltemp))
+ parameterline=13;
+ else{
+ parameterline=14;
}
- if (feof(ficpar))
+ switch (parameterline){
+ case 11:
+ if((num_filled=sscanf(line,"prevforecast=%d starting-proj-date=%lf/%lf/%lf final-proj-date=%lf/%lf/%lf mobil_average=%d\n",&prevfcast,&jproj1,&mproj1,&anproj1,&jproj2,&mproj2,&anproj2,&mobilavproj)) !=EOF){
+ if (num_filled != 8) {
+ printf("Error: Not 8 (data)parameters in line but %d, for example:prevforecast=1 starting-proj-date=1/1/1990 final-proj-date=1/1/2000 mobil_average=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ fprintf(ficlog,"Error: Not 8 (data)parameters in line but %d, for example:prevforecast=1 starting-proj-date=1/1/1990 final-proj-date=1/1/2000 mov_average=0\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ goto end;
+ }
+ fprintf(ficparo,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
+ printf("prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
+ fprintf(ficlog,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
+ fprintf(ficres,"prevforecast=%d starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mobil_average=%d\n",prevfcast,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilavproj);
+ /* day and month of proj2 are not used but only year anproj2.*/
+ }
break;
- else{ /* Processess output results for this combination of covariate values */
- }
- } /* end while */
-
-
+ case 12:
+ /*fscanf(ficpar,"backcast=%d starting-back-date=%lf/%lf/%lf final-back-date=%lf/%lf/%lf mobil_average=%d\n",&backcast,&jback1,&mback1,&anback1,&jback2,&mback2,&anback2,&mobilavproj);*/
+ if((num_filled=sscanf(line,"backcast=%d starting-back-date=%lf/%lf/%lf final-back-date=%lf/%lf/%lf mobil_average=%d\n",&backcast,&jback1,&mback1,&anback1,&jback2,&mback2,&anback2,&mobilavproj)) !=EOF){
+ if (num_filled != 8) {
+ printf("Error: Not 8 (data)parameters in line but %d, for example:backcast=1 starting-back-date=1/1/1990 finloal-back-date=1/1/1970 mobil_average=1\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ fprintf(ficlog,"Error: Not 8 (data)parameters in line but %d, for example:backcast=1 starting-back-date=1/1/1990 finloal-back-date=1/1/1970 mobil_average=1\n, your line=%s . Probably you are running an older format.\n",num_filled,line);
+ goto end;
+ }
+ printf("backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ fprintf(ficparo,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ fprintf(ficlog,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ fprintf(ficres,"backcast=%d starting-back-date=%.lf/%.lf/%.lf final-back-date=%.lf/%.lf/%.lf mobil_average=%d\n",backcast,jback1,mback1,anback1,jback2,mback2,anback2,mobilavproj);
+ /* day and month of proj2 are not used but only year anproj2.*/
+ }
+ break;
+ case 13:
+ if((num_filled=sscanf(line,"result:%[^\n]\n",resultline)) !=EOF){
+ if (num_filled == 0){
+ resultline[0]='\0';
+ printf("Warning %d: no result line! It should be at minimum 'result: V2=0 V1=1 or result:.\n%s\n", num_filled, line);
+ fprintf(ficlog,"Warning %d: no result line! It should be at minimum 'result: V2=0 V1=1 or result:.\n%s\n", num_filled, line);
+ break;
+ } else if (num_filled != 1){
+ printf("ERROR %d: result line! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",num_filled, line);
+ fprintf(ficlog,"ERROR %d: result line! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",num_filled, line);
+ }
+ nresult++; /* Sum of resultlines */
+ printf("Result %d: result=%s\n",nresult, resultline);
+ if(nresult > MAXRESULTLINES){
+ printf("ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
+ fprintf(ficlog,"ERROR: Current version of IMaCh limits the number of resultlines to %d, you used %d\n",MAXRESULTLINES,nresult);
+ goto end;
+ }
+ decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
+ fprintf(ficparo,"result: %s\n",resultline);
+ fprintf(ficres,"result: %s\n",resultline);
+ fprintf(ficlog,"result: %s\n",resultline);
+ break;
+ case 14:
+ if(ncovmodel >2 && nresult==0 ){
+ printf("ERROR: no result lines! It should be at minimum 'result: V2=0 V1=1 or result:.' %s\n",line);
+ goto end;
+ }
+ break;
+ default:
+ nresult=1;
+ decoderesult(".",nresult ); /* No covariate */
+ }
+ } /* End switch parameterline */
+ }while(endishere==0); /* End do */
/* freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint); */
/* ,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2); */
@@ -11127,7 +11298,7 @@ Please run with mle=-1 to get a correct
printinggnuplot(fileresu, optionfilefiname,ageminpar,agemaxpar,fage, prevfcast, backcast, pathc,p);
}
printinghtml(fileresu,title,datafile, firstpass, lastpass, stepm, weightopt, \
- model,imx,jmin,jmax,jmean,rfileres,popforecast,prevfcast,backcast, estepm, \
+ model,imx,jmin,jmax,jmean,rfileres,popforecast,mobilav,prevfcast,mobilavproj,backcast, estepm, \
jprev1,mprev1,anprev1,dateprev1,jprev2,mprev2,anprev2,dateprev2);
/*------------ free_vector -------------*/
@@ -11180,6 +11351,7 @@ Please run with mle=-1 to get a correct
mobaverage=mobaverages;
if (mobilav!=0) {
printf("Movingaveraging observed prevalence\n");
+ fprintf(ficlog,"Movingaveraging observed prevalence\n");
if (movingaverage(probs, ageminpar, agemaxpar, mobaverage, mobilav)!=0){
fprintf(ficlog," Error in movingaverage mobilav=%d\n",mobilav);
printf(" Error in movingaverage mobilav=%d\n",mobilav);
@@ -11189,6 +11361,7 @@ Please run with mle=-1 to get a correct
/* prevalence(probs, ageminpar, agemaxpar, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass); */
else if (mobilavproj !=0) {
printf("Movingaveraging projected observed prevalence\n");
+ fprintf(ficlog,"Movingaveraging projected observed prevalence\n");
if (movingaverage(probs, ageminpar, agemaxpar, mobaverage, mobilavproj)!=0){
fprintf(ficlog," Error in movingaverage mobilavproj=%d\n",mobilavproj);
printf(" Error in movingaverage mobilavproj=%d\n",mobilavproj);
@@ -11251,7 +11424,7 @@ Please run with mle=-1 to get a correct
for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
fprintf(ficreseij,"\n#****** ");
printf("\n#****** ");
@@ -11324,7 +11497,7 @@ Please run with mle=-1 to get a correct
for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
printf("\n#****** Result for:");
fprintf(ficrest,"\n#****** Result for:");
@@ -11465,7 +11638,7 @@ Please run with mle=-1 to get a correct
for(nres=1; nres <= nresult; nres++) /* For each resultline */
for(k=1; k<=i1;k++){
- if(TKresult[nres]!= k)
+ if(i1 != 1 && TKresult[nres]!= k)
continue;
fprintf(ficresvpl,"\n#****** ");
printf("\n#****** ");
![%s_%d-2-%d.svg](\"%s_%d-2-%d.svg\"){kind=link}