|
|
| version 1.283, 2018/04/19 14:49:16 | version 1.293, 2019/05/09 15:17:34 |
|---|---|
| Line 1 | Line 1 |
| /* $Id$ | /* $Id$ |
| $State$ | $State$ |
| $Log$ | $Log$ |
| Revision 1.293 2019/05/09 15:17:34 brouard | |
| *** empty log message *** | |
| Revision 1.292 2019/05/09 14:17:20 brouard | |
| Summary: Some updates | |
| Revision 1.291 2019/05/09 13:44:18 brouard | |
| Summary: Before ncovmax | |
| Revision 1.290 2019/05/09 13:39:37 brouard | |
| Summary: 0.99r18 unlimited number of individuals | |
| The number n which was limited to 20,000 cases is now unlimited, from firstobs to lastobs. If the number is too for the virtual memory, probably an error will occur. | |
| Revision 1.289 2018/12/13 09:16:26 brouard | |
| Summary: Bug for young ages (<-30) will be in r17 | |
| Revision 1.288 2018/05/02 20:58:27 brouard | |
| Summary: Some bugs fixed | |
| Revision 1.287 2018/05/01 17:57:25 brouard | |
| Summary: Bug fixed by providing frequencies only for non missing covariates | |
| Revision 1.286 2018/04/27 14:27:04 brouard | |
| Summary: some minor bugs | |
| Revision 1.285 2018/04/21 21:02:16 brouard | |
| Summary: Some bugs fixed, valgrind tested | |
| Revision 1.284 2018/04/20 05:22:13 brouard | |
| Summary: Computing mean and stdeviation of fixed quantitative variables | |
| Revision 1.283 2018/04/19 14:49:16 brouard | Revision 1.283 2018/04/19 14:49:16 brouard |
| Summary: Some minor bugs fixed | Summary: Some minor bugs fixed |
| Line 1031 typedef struct { | Line 1063 typedef struct { |
| #define NINTERVMAX 8 | #define NINTERVMAX 8 |
| #define NLSTATEMAX 8 /**< Maximum number of live states (for func) */ | #define NLSTATEMAX 8 /**< Maximum number of live states (for func) */ |
| #define NDEATHMAX 8 /**< Maximum number of dead states (for func) */ | #define NDEATHMAX 8 /**< Maximum number of dead states (for func) */ |
| #define NCOVMAX 20 /**< Maximum number of covariates, including generated covariates V1*V2 */ | #define NCOVMAX 20 /**< Maximum number of covariates, including generated covariates V1*V2 */ |
| #define codtabm(h,k) (1 & (h-1) >> (k-1))+1 | #define codtabm(h,k) (1 & (h-1) >> (k-1))+1 |
| /*#define decodtabm(h,k,cptcoveff)= (h <= (1<<cptcoveff)?(((h-1) >> (k-1)) & 1) +1 : -1)*/ | /*#define decodtabm(h,k,cptcoveff)= (h <= (1<<cptcoveff)?(((h-1) >> (k-1)) & 1) +1 : -1)*/ |
| #define decodtabm(h,k,cptcoveff) (((h-1) >> (k-1)) & 1) +1 | #define decodtabm(h,k,cptcoveff) (((h-1) >> (k-1)) & 1) +1 |
| #define MAXN 20000 | /*#define MAXN 20000 */ /* Should by replaced by nobs, real number of observations and unlimited */ |
| #define YEARM 12. /**< Number of months per year */ | #define YEARM 12. /**< Number of months per year */ |
| /* #define AGESUP 130 */ | /* #define AGESUP 130 */ |
| #define AGESUP 150 | /* #define AGESUP 150 */ |
| #define AGESUP 200 | |
| #define AGEINF 0 | #define AGEINF 0 |
| #define AGEMARGE 25 /* Marge for agemin and agemax for(iage=agemin-AGEMARGE; iage <= agemax+3+AGEMARGE; iage++) */ | #define AGEMARGE 25 /* Marge for agemin and agemax for(iage=agemin-AGEMARGE; iage <= agemax+3+AGEMARGE; iage++) */ |
| #define AGEBASE 40 | #define AGEBASE 40 |
| Line 1082 int nqfveff=0; /**< nqfveff Number of Qu | Line 1115 int nqfveff=0; /**< nqfveff Number of Qu |
| int ntveff=0; /**< ntveff number of effective time varying variables */ | int ntveff=0; /**< ntveff number of effective time varying variables */ |
| int nqtveff=0; /**< ntqveff number of effective time varying quantitative variables */ | int nqtveff=0; /**< ntqveff number of effective time varying quantitative variables */ |
| int cptcov=0; /* Working variable */ | int cptcov=0; /* Working variable */ |
| int nobs=10; /* Number of observations in the data lastobs-firstobs */ | |
| int ncovcombmax=NCOVMAX; /* Maximum calculated number of covariate combination = pow(2, cptcoveff) */ | int ncovcombmax=NCOVMAX; /* Maximum calculated number of covariate combination = pow(2, cptcoveff) */ |
| int npar=NPARMAX; | int npar=NPARMAX; |
| int nlstate=2; /* Number of live states */ | int nlstate=2; /* Number of live states */ |
| Line 2565 void powell(double p[], double **xi, int | Line 2599 void powell(double p[], double **xi, int |
| double **newm; | double **newm; |
| double agefin, delaymax=200. ; /* 100 Max number of years to converge */ | double agefin, delaymax=200. ; /* 100 Max number of years to converge */ |
| int ncvloop=0; | int ncvloop=0; |
| int first=0; | |
| min=vector(1,nlstate); | min=vector(1,nlstate); |
| max=vector(1,nlstate); | max=vector(1,nlstate); |
| Line 2661 void powell(double p[], double **xi, int | Line 2696 void powell(double p[], double **xi, int |
| free_vector(meandiff,1,nlstate); | free_vector(meandiff,1,nlstate); |
| return prlim; | return prlim; |
| } | } |
| } /* age loop */ | } /* agefin loop */ |
| /* After some age loop it doesn't converge */ | /* After some age loop it doesn't converge */ |
| printf("Warning: the stable prevalence at age %d did not converge with the required precision (%g > ftolpl=%g) within %.0f years. Try to lower 'ftolpl'. \n\ | if(!first){ |
| Earliest age to start was %d-%d=%d, ncvloop=%d, ncvyear=%d\n", (int)age, maxmax, ftolpl, delaymax, (int)age, (int)delaymax, (int)agefin, ncvloop, *ncvyear); | first=1; |
| printf("Warning: the stable prevalence at age %d did not converge with the required precision (%g > ftolpl=%g) within %.d years and %d loops. Try to lower 'ftolpl'. Youngest age to start was %d=(%d-%d). Others in log file only...\n", (int)age, maxmax, ftolpl, *ncvyear, ncvloop, (int)(agefin+stepm/YEARM), (int)(age-stepm/YEARM), (int)delaymax); | |
| } | |
| fprintf(ficlog, "Warning: the stable prevalence at age %d did not converge with the required precision (%g > ftolpl=%g) within %.d years and %d loops. Try to lower 'ftolpl'. Youngest age to start was %d=(%d-%d).\n", (int)age, maxmax, ftolpl, *ncvyear, ncvloop, (int)(agefin+stepm/YEARM), (int)(age-stepm/YEARM), (int)delaymax); | |
| /* Try to lower 'ftol', for example from 1.e-8 to 6.e-9.\n", ftolpl, (int)age, (int)delaymax, (int)agefin, ncvloop, (int)age-(int)agefin); */ | /* Try to lower 'ftol', for example from 1.e-8 to 6.e-9.\n", ftolpl, (int)age, (int)delaymax, (int)agefin, ncvloop, (int)age-(int)agefin); */ |
| free_vector(min,1,nlstate); | free_vector(min,1,nlstate); |
| free_vector(max,1,nlstate); | free_vector(max,1,nlstate); |
| Line 2730 Earliest age to start was %d-%d=%d, ncvl | Line 2769 Earliest age to start was %d-%d=%d, ncvl |
| /* Even if hstepm = 1, at least one multiplication by the unit matrix */ | /* Even if hstepm = 1, at least one multiplication by the unit matrix */ |
| /* Start at agefin= age, computes the matrix of passage and loops decreasing agefin until convergence is reached */ | /* Start at agefin= age, computes the matrix of passage and loops decreasing agefin until convergence is reached */ |
| /* for(agefin=age+stepm/YEARM; agefin<=age+delaymax; agefin=agefin+stepm/YEARM){ /\* A changer en age *\/ */ | /* for(agefin=age+stepm/YEARM; agefin<=age+delaymax; agefin=agefin+stepm/YEARM){ /\* A changer en age *\/ */ |
| for(agefin=age; agefin<AGESUP; agefin=agefin+stepm/YEARM){ /* A changer en age */ | /* for(agefin=age; agefin<AGESUP; agefin=agefin+stepm/YEARM){ /\* A changer en age *\/ */ |
| for(agefin=age; agefin<FMIN(AGESUP,age+delaymax); agefin=agefin+stepm/YEARM){ /* A changer en age */ | |
| ncvloop++; | ncvloop++; |
| newm=savm; /* oldm should be kept from previous iteration or unity at start */ | newm=savm; /* oldm should be kept from previous iteration or unity at start */ |
| /* newm points to the allocated table savm passed by the function it can be written, savm could be reallocated */ | /* newm points to the allocated table savm passed by the function it can be written, savm could be reallocated */ |
| Line 2844 Earliest age to start was %d-%d=%d, ncvl | Line 2884 Earliest age to start was %d-%d=%d, ncvl |
| free_vector(meandiff,1,nlstate); | free_vector(meandiff,1,nlstate); |
| return bprlim; | return bprlim; |
| } | } |
| } /* age loop */ | } /* agefin loop */ |
| /* After some age loop it doesn't converge */ | /* After some age loop it doesn't converge */ |
| if(first){ | if(!first){ |
| first=1; | first=1; |
| printf("Warning: the back stable prevalence at age %d did not converge with the required precision (%g > ftolpl=%g) within %.0f years. Try to lower 'ftolpl'. Others in log file only...\n\ | printf("Warning: the back stable prevalence at age %d did not converge with the required precision (%g > ftolpl=%g) within %.0f years. Try to lower 'ftolpl'. Others in log file only...\n\ |
| Oldest age to start was %d-%d=%d, ncvloop=%d, ncvyear=%d\n", (int)age, maxmax, ftolpl, delaymax, (int)age, (int)delaymax, (int)agefin, ncvloop, *ncvyear); | Oldest age to start was %d-%d=%d, ncvloop=%d, ncvyear=%d\n", (int)age, maxmax, ftolpl, delaymax, (int)age, (int)delaymax, (int)agefin, ncvloop, *ncvyear); |
| Line 2929 double **pmij(double **ps, double *cov, | Line 2969 double **pmij(double **ps, double *cov, |
| ps[ii][ii]=1; | ps[ii][ii]=1; |
| } | } |
| } | } |
| /* Added for backcast */ /* Transposed matrix too */ | |
| for(jj=1; jj<= nlstate+ndeath; jj++){ | |
| s1=0.; | |
| for(ii=1; ii<= nlstate+ndeath; ii++){ | |
| s1+=ps[ii][jj]; | |
| } | |
| for(ii=1; ii<= nlstate; ii++){ | |
| ps[ii][jj]=ps[ii][jj]/s1; | |
| } | |
| } | |
| /* Transposition */ | |
| for(jj=1; jj<= nlstate+ndeath; jj++){ | |
| for(ii=jj; ii<= nlstate+ndeath; ii++){ | |
| s1=ps[ii][jj]; | |
| ps[ii][jj]=ps[jj][ii]; | |
| ps[jj][ii]=s1; | |
| } | |
| } | |
| /* for(ii=1; ii<= nlstate+ndeath; ii++){ */ | /* for(ii=1; ii<= nlstate+ndeath; ii++){ */ |
| /* for(jj=1; jj<= nlstate+ndeath; jj++){ */ | /* for(jj=1; jj<= nlstate+ndeath; jj++){ */ |
| /* printf(" pmij ps[%d][%d]=%lf ",ii,jj,ps[ii][jj]); */ | /* printf(" pmij ps[%d][%d]=%lf ",ii,jj,ps[ii][jj]); */ |
| Line 2958 double **pmij(double **ps, double *cov, | Line 3014 double **pmij(double **ps, double *cov, |
| double **out, **pmij(); | double **out, **pmij(); |
| double sumnew=0.; | double sumnew=0.; |
| double agefin; | double agefin; |
| double k3=0.; /* constant of the w_x diagonal matrixe (in order for B to sum to 1 even for death state) */ | double k3=0.; /* constant of the w_x diagonal matrix (in order for B to sum to 1 even for death state) */ |
| double **dnewm, **dsavm, **doldm; | double **dnewm, **dsavm, **doldm; |
| double **bbmij; | double **bbmij; |
| Line 2977 double **pmij(double **ps, double *cov, | Line 3033 double **pmij(double **ps, double *cov, |
| /* outputs pmmij which is a stochastic matrix in row */ | /* outputs pmmij which is a stochastic matrix in row */ |
| /* Diag(w_x) */ | /* Diag(w_x) */ |
| /* Problem with prevacurrent which can be zero */ | /* Rescaling the cross-sectional prevalence: Problem with prevacurrent which can be zero */ |
| sumnew=0.; | sumnew=0.; |
| /*for (ii=1;ii<=nlstate+ndeath;ii++){*/ | /*for (ii=1;ii<=nlstate+ndeath;ii++){*/ |
| for (ii=1;ii<=nlstate;ii++){ /* Only on live states */ | for (ii=1;ii<=nlstate;ii++){ /* Only on live states */ |
| Line 3004 double **pmij(double **ps, double *cov, | Line 3060 double **pmij(double **ps, double *cov, |
| } | } |
| /* End doldm, At the end doldm is diag[(w_i)] */ | /* End doldm, At the end doldm is diag[(w_i)] */ |
| /* left Product of this diag matrix by pmmij=Px (dnewm=dsavm*doldm) */ | /* Left product of this diag matrix by pmmij=Px (dnewm=dsavm*doldm): diag[(w_i)*Px */ |
| bbmij=matprod2(dnewm, doldm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, pmmij); /* Bug Valgrind */ | bbmij=matprod2(dnewm, doldm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, pmmij); /* was a Bug Valgrind */ |
| /* Diag(Sum_i w^i_x p^ij_x */ | /* Diag(Sum_i w^i_x p^ij_x, should be the prevalence at age x+stepm */ |
| /* w1 p11 + w2 p21 only on live states N1./N..*N11/N1. + N2./N..*N21/N2.=(N11+N21)/N..=N.1/N.. */ | /* w1 p11 + w2 p21 only on live states N1./N..*N11/N1. + N2./N..*N21/N2.=(N11+N21)/N..=N.1/N.. */ |
| for (j=1;j<=nlstate+ndeath;j++){ | for (j=1;j<=nlstate+ndeath;j++){ |
| sumnew=0.; | sumnew=0.; |
| Line 3025 double **pmij(double **ps, double *cov, | Line 3081 double **pmij(double **ps, double *cov, |
| } /*End ii */ | } /*End ii */ |
| } /* End j, At the end dsavm is diag[1/(w_1p1i+w_2 p2i)] for ALL states even if the sum is only for live states */ | } /* End j, At the end dsavm is diag[1/(w_1p1i+w_2 p2i)] for ALL states even if the sum is only for live states */ |
| ps=matprod2(ps, dnewm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, dsavm); /* Bug Valgrind */ | ps=matprod2(ps, dnewm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, dsavm); /* was a Bug Valgrind */ |
| /* ps is now diag[w_i] * Px * diag [1/(w_1p1i+w_2 p2i)] */ | /* ps is now diag[w_i] * Px * diag [1/(w_1p1i+w_2 p2i)] */ |
| /* end bmij */ | /* end bmij */ |
| return ps; /*pointer is unchanged */ | return ps; /*pointer is unchanged */ |
| Line 3857 return -l; | Line 3913 return -l; |
| /*************** function likelione ***********/ | /*************** function likelione ***********/ |
| void likelione(FILE *ficres,double p[], int npar, int nlstate, int *globpri, long *ipmx, double *sw, double *fretone, double (*funcone)(double [])) | void likelione(FILE *ficres,double p[], int npar, int nlstate, int *globpri, long *ipmx, double *sw, double *fretone, double (*func)(double [])) |
| { | { |
| /* This routine should help understanding what is done with | /* This routine should help understanding what is done with |
| the selection of individuals/waves and | the selection of individuals/waves and |
| Line 3881 void likelione(FILE *ficres,double p[], | Line 3937 void likelione(FILE *ficres,double p[], |
| fprintf(ficresilk," -2*gipw/gsw*weight*ll(total)\n"); | fprintf(ficresilk," -2*gipw/gsw*weight*ll(total)\n"); |
| } | } |
| *fretone=(*funcone)(p); | *fretone=(*func)(p); |
| if(*globpri !=0){ | if(*globpri !=0){ |
| fclose(ficresilk); | fclose(ficresilk); |
| if (mle ==0) | if (mle ==0) |
| Line 4371 void freqsummary(char fileres[], double | Line 4427 void freqsummary(char fileres[], double |
| double ***freq; /* Frequencies */ | double ***freq; /* Frequencies */ |
| double *x, *y, a=0.,b=0.,r=1., sa=0., sb=0.; /* for regression, y=b+m*x and r is the correlation coefficient */ | double *x, *y, a=0.,b=0.,r=1., sa=0., sb=0.; /* for regression, y=b+m*x and r is the correlation coefficient */ |
| int no=0, linreg(int ifi, int ila, int *no, const double x[], const double y[], double* a, double* b, double* r, double* sa, double * sb); | int no=0, linreg(int ifi, int ila, int *no, const double x[], const double y[], double* a, double* b, double* r, double* sa, double * sb); |
| double *meanq, *idq; | double *meanq, *stdq, *idq; |
| double **meanqt; | double **meanqt; |
| double *pp, **prop, *posprop, *pospropt; | double *pp, **prop, *posprop, *pospropt; |
| double pos=0., posproptt=0., pospropta=0., k2, dateintsum=0,k2cpt=0; | double pos=0., posproptt=0., pospropta=0., k2, dateintsum=0,k2cpt=0; |
| Line 4384 void freqsummary(char fileres[], double | Line 4440 void freqsummary(char fileres[], double |
| pospropt=vector(1,nlstate); /* Counting the number of transition starting from a live state */ | pospropt=vector(1,nlstate); /* Counting the number of transition starting from a live state */ |
| /* prop=matrix(1,nlstate,iagemin,iagemax+3); */ | /* prop=matrix(1,nlstate,iagemin,iagemax+3); */ |
| meanq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */ | meanq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */ |
| stdq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */ | |
| idq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */ | idq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */ |
| meanqt=matrix(1,lastpass,1,nqtveff); | meanqt=matrix(1,lastpass,1,nqtveff); |
| strcpy(fileresp,"P_"); | strcpy(fileresp,"P_"); |
| Line 4494 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4551 Title=%s <br>Datafile=%s Firstpass=%d La |
| pospropt[i]=0; | pospropt[i]=0; |
| } | } |
| for (z1=1; z1<= nqfveff; z1++) { /* zeroing for each combination j1 as well as for the total */ | for (z1=1; z1<= nqfveff; z1++) { /* zeroing for each combination j1 as well as for the total */ |
| idq[z1]+=0.; | idq[z1]=0.; |
| meanq[z1]+=0.; | meanq[z1]=0.; |
| stdq[z1]=0.; | |
| } | } |
| /* for (z1=1; z1<= nqtveff; z1++) { */ | /* for (z1=1; z1<= nqtveff; z1++) { */ |
| /* for(m=1;m<=lastpass;m++){ */ | /* for(m=1;m<=lastpass;m++){ */ |
| Line 4511 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4569 Title=%s <br>Datafile=%s Firstpass=%d La |
| if(j !=0){ | if(j !=0){ |
| if(anyvaryingduminmodel==0){ /* If All fixed covariates */ | if(anyvaryingduminmodel==0){ /* If All fixed covariates */ |
| if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */ | if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */ |
| /* for (z1=1; z1<= nqfveff; z1++) { */ | |
| /* meanq[z1]+=coqvar[Tvar[z1]][iind]; /\* Computes mean of quantitative with selected filter *\/ */ | |
| /* } */ | |
| for (z1=1; z1<=cptcoveff; z1++) { /* loops on covariates in the model */ | for (z1=1; z1<=cptcoveff; z1++) { /* loops on covariates in the model */ |
| /* if(Tvaraff[z1] ==-20){ */ | /* if(Tvaraff[z1] ==-20){ */ |
| /* /\* sumnew+=cotvar[mw[mi][iind]][z1][iind]; *\/ */ | /* /\* sumnew+=cotvar[mw[mi][iind]][z1][iind]; *\/ */ |
| Line 4534 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4589 Title=%s <br>Datafile=%s Firstpass=%d La |
| }/* end j==0 */ | }/* end j==0 */ |
| if (bool==1){ /* We selected an individual iind satisfying combination j1 (V4=1 V3=0) or all fixed covariates */ | if (bool==1){ /* We selected an individual iind satisfying combination j1 (V4=1 V3=0) or all fixed covariates */ |
| /* for(m=firstpass; m<=lastpass; m++){ */ | /* for(m=firstpass; m<=lastpass; m++){ */ |
| for(mi=1; mi<wav[iind];mi++){ /* For that wave */ | for(mi=1; mi<wav[iind];mi++){ /* For each wave */ |
| m=mw[mi][iind]; | m=mw[mi][iind]; |
| if(j!=0){ | if(j!=0){ |
| if(anyvaryingduminmodel==1){ /* Some are varying covariates */ | if(anyvaryingduminmodel==1){ /* Some are varying covariates */ |
| Line 4554 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4609 Title=%s <br>Datafile=%s Firstpass=%d La |
| }/* Some are varying covariates, we tried to speed up if all fixed covariates in the model, avoiding waves loop */ | }/* Some are varying covariates, we tried to speed up if all fixed covariates in the model, avoiding waves loop */ |
| } /* end j==0 */ | } /* end j==0 */ |
| /* bool =0 we keep that guy which corresponds to the combination of dummy values */ | /* bool =0 we keep that guy which corresponds to the combination of dummy values */ |
| if(bool==1){ | if(bool==1){ /*Selected */ |
| /* dh[m][iind] or dh[mw[mi][iind]][iind] is the delay between two effective (mi) waves m=mw[mi][iind] | /* dh[m][iind] or dh[mw[mi][iind]][iind] is the delay between two effective (mi) waves m=mw[mi][iind] |
| and mw[mi+1][iind]. dh depends on stepm. */ | and mw[mi+1][iind]. dh depends on stepm. */ |
| agebegin=agev[m][iind]; /* Age at beginning of wave before transition*/ | agebegin=agev[m][iind]; /* Age at beginning of wave before transition*/ |
| Line 4572 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4627 Title=%s <br>Datafile=%s Firstpass=%d La |
| if(s[m][iind]==-1) | if(s[m][iind]==-1) |
| printf(" num=%ld m=%d, iind=%d s1=%d s2=%d agev at m=%d agebegin=%.2f ageend=%.2f, agemed=%d\n", num[iind], m, iind,s[m][iind],s[m+1][iind], (int)agev[m][iind],agebegin, ageend, (int)((agebegin+ageend)/2.)); | printf(" num=%ld m=%d, iind=%d s1=%d s2=%d agev at m=%d agebegin=%.2f ageend=%.2f, agemed=%d\n", num[iind], m, iind,s[m][iind],s[m+1][iind], (int)agev[m][iind],agebegin, ageend, (int)((agebegin+ageend)/2.)); |
| freq[s[m][iind]][s[m+1][iind]][(int)agev[m][iind]] += weight[iind]; /* At age of beginning of transition, where status is known */ | freq[s[m][iind]][s[m+1][iind]][(int)agev[m][iind]] += weight[iind]; /* At age of beginning of transition, where status is known */ |
| for (z1=1; z1<= nqfveff; z1++) { /* Quantitative variables, calculating mean */ | |
| idq[z1]=idq[z1]+weight[iind]; | |
| meanq[z1]+=covar[ncovcol+z1][iind]*weight[iind]; /* Computes mean of quantitative with selected filter */ | |
| stdq[z1]+=covar[ncovcol+z1][iind]*covar[ncovcol+z1][iind]*weight[iind]*weight[iind]; /* *weight[iind];*/ /* Computes mean of quantitative with selected filter */ | |
| } | |
| /* if((int)agev[m][iind] == 55) */ | /* if((int)agev[m][iind] == 55) */ |
| /* printf("j=%d, j1=%d Age %d, iind=%d, num=%09ld m=%d\n",j,j1,(int)agev[m][iind],iind, num[iind],m); */ | /* printf("j=%d, j1=%d Age %d, iind=%d, num=%09ld m=%d\n",j,j1,(int)agev[m][iind],iind, num[iind],m); */ |
| /* freq[s[m][iind]][s[m+1][iind]][(int)((agebegin+ageend)/2.)] += weight[iind]; */ | /* freq[s[m][iind]][s[m+1][iind]][(int)((agebegin+ageend)/2.)] += weight[iind]; */ |
| freq[s[m][iind]][s[m+1][iind]][iagemax+3] += weight[iind]; /* Total is in iagemax+3 *//* At age of beginning of transition, where status is known */ | freq[s[m][iind]][s[m+1][iind]][iagemax+3] += weight[iind]; /* Total is in iagemax+3 *//* At age of beginning of transition, where status is known */ |
| } | } |
| for (z1=1; z1<= nqfveff; z1++) { | |
| idq[z1]++; | |
| meanq[z1]+=covar[ncovcol+z1][iind]; /* *weight[iind];*/ /* Computes mean of quantitative with selected filter */ | |
| /* meanq[z1]+=coqvar[Tvar[z1]][iind]; /\* Computes mean of quantitative with selected filter *\/ */ | |
| } | |
| } /* end if between passes */ | } /* end if between passes */ |
| if ((agev[m][iind]>1) && (agev[m][iind]< (iagemax+3)) && (anint[m][iind]!=9999) && (mint[m][iind]!=99) && (j==0)) { | if ((agev[m][iind]>1) && (agev[m][iind]< (iagemax+3)) && (anint[m][iind]!=9999) && (mint[m][iind]!=99) && (j==0)) { |
| dateintsum=dateintsum+k2; /* on all covariates ?*/ | dateintsum=dateintsum+k2; /* on all covariates ?*/ |
| Line 4592 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4647 Title=%s <br>Datafile=%s Firstpass=%d La |
| bool=1; | bool=1; |
| }/* end bool 2 */ | }/* end bool 2 */ |
| } /* end m */ | } /* end m */ |
| /* for (z1=1; z1<= nqfveff; z1++) { /\* Quantitative variables, calculating mean *\/ */ | |
| /* idq[z1]=idq[z1]+weight[iind]; */ | |
| /* meanq[z1]+=covar[ncovcol+z1][iind]*weight[iind]; /\* Computes mean of quantitative with selected filter *\/ */ | |
| /* stdq[z1]+=covar[ncovcol+z1][iind]*covar[ncovcol+z1][iind]*weight[iind]*weight[iind]; /\* *weight[iind];*\/ /\* Computes mean of quantitative with selected filter *\/ */ | |
| /* } */ | |
| } /* end bool */ | } /* end bool */ |
| } /* end iind = 1 to imx */ | } /* end iind = 1 to imx */ |
| /* prop[s][age] is feeded for any initial and valid live state as well as | /* prop[s][age] is feeded for any initial and valid live state as well as |
| Line 4629 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4689 Title=%s <br>Datafile=%s Firstpass=%d La |
| fprintf(ficresphtmfr, "**********</h3>\n"); | fprintf(ficresphtmfr, "**********</h3>\n"); |
| fprintf(ficlog, "**********\n"); | fprintf(ficlog, "**********\n"); |
| } | } |
| /* | /* |
| Printing means of quantitative variables if any | Printing means of quantitative variables if any |
| */ | */ |
| for (z1=1; z1<= nqfveff; z1++) { | for (z1=1; z1<= nqfveff; z1++) { |
| fprintf(ficresphtmfr,"V quantitative id %d, number of idividuals= %f, sum=%f", z1, idq[z1], meanq[z1]); | fprintf(ficlog,"Mean of fixed quantitative variable V%d on %.0f individuals sum=%f", ncovcol+z1, idq[z1], meanq[z1]); |
| fprintf(ficresphtmfr,", mean=%f<p>\n",meanq[z1]/idq[z1]); | fprintf(ficlog,", mean=%.3g\n",meanq[z1]/idq[z1]); |
| } | if(weightopt==1){ |
| /* for (z1=1; z1<= nqtveff; z1++) { */ | printf(" Weighted mean and standard deviation of"); |
| /* for(m=1;m<=lastpass;m++){ */ | fprintf(ficlog," Weighted mean and standard deviation of"); |
| /* fprintf(ficresphtmfr,"V quantitative id %d, pass id=%d, mean=%f<p>\n", z1, m, meanqt[m][z1]); */ | fprintf(ficresphtmfr," Weighted mean and standard deviation of"); |
| /* } */ | } |
| /* } */ | printf(" fixed quantitative variable V%d on %.0f representatives of the population : %6.3g (%6.3g)\n", ncovcol+z1, idq[z1],meanq[z1]/idq[z1], sqrt((stdq[z1]-meanq[z1]*meanq[z1]/idq[z1])/idq[z1])); |
| fprintf(ficlog," fixed quantitative variable V%d on %.0f representatives of the population : %6.3g (%6.3g)\n", ncovcol+z1, idq[z1],meanq[z1]/idq[z1], sqrt((stdq[z1]-meanq[z1]*meanq[z1]/idq[z1])/idq[z1])); | |
| fprintf(ficresphtmfr," fixed quantitative variable V%d on %.0f representatives of the population : %6.3g (%6.3g)<p>\n", ncovcol+z1, idq[z1],meanq[z1]/idq[z1], sqrt((stdq[z1]-meanq[z1]*meanq[z1]/idq[z1])/idq[z1])); | |
| } | |
| /* for (z1=1; z1<= nqtveff; z1++) { */ | |
| /* for(m=1;m<=lastpass;m++){ */ | |
| /* fprintf(ficresphtmfr,"V quantitative id %d, pass id=%d, mean=%f<p>\n", z1, m, meanqt[m][z1]); */ | |
| /* } */ | |
| /* } */ | |
| fprintf(ficresphtm,"<table style=\"text-align:center; border: 1px solid\">"); | fprintf(ficresphtm,"<table style=\"text-align:center; border: 1px solid\">"); |
| if((cptcoveff==0 && nj==1)|| nj==2 ) /* no covariate and first pass */ | if((cptcoveff==0 && nj==1)|| nj==2 ) /* no covariate and first pass */ |
| Line 4876 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4944 Title=%s <br>Datafile=%s Firstpass=%d La |
| fprintf(ficlog,"\n"); | fprintf(ficlog,"\n"); |
| } | } |
| } | } |
| } | } /* end of state i */ |
| printf("#Freqsummary\n"); | printf("#Freqsummary\n"); |
| fprintf(ficlog,"\n"); | fprintf(ficlog,"\n"); |
| for(s1=-1; s1 <=nlstate+ndeath; s1++){ | for(s1=-1; s1 <=nlstate+ndeath; s1++){ |
| Line 4924 Title=%s <br>Datafile=%s Firstpass=%d La | Line 4992 Title=%s <br>Datafile=%s Firstpass=%d La |
| fclose(ficresphtmfr); | fclose(ficresphtmfr); |
| free_vector(idq,1,nqfveff); | free_vector(idq,1,nqfveff); |
| free_vector(meanq,1,nqfveff); | free_vector(meanq,1,nqfveff); |
| free_vector(stdq,1,nqfveff); | |
| free_matrix(meanqt,1,lastpass,1,nqtveff); | free_matrix(meanqt,1,lastpass,1,nqtveff); |
| free_vector(x, iagemin-AGEMARGE, iagemax+4+AGEMARGE); | free_vector(x, iagemin-AGEMARGE, iagemax+4+AGEMARGE); |
| free_vector(y, iagemin-AGEMARGE, iagemax+4+AGEMARGE); | free_vector(y, iagemin-AGEMARGE, iagemax+4+AGEMARGE); |
| Line 5030 void prevalence(double ***probs, double | Line 5099 void prevalence(double ***probs, double |
| /*j=cptcoveff;*/ | /*j=cptcoveff;*/ |
| if (cptcovn<1) {j=1;ncodemax[1]=1;} | if (cptcovn<1) {j=1;ncodemax[1]=1;} |
| first=1; | first=0; |
| for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of covariate */ | for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of covariate */ |
| for (i=1; i<=nlstate; i++) | for (i=1; i<=nlstate; i++) |
| for(iage=iagemin-AGEMARGE; iage <= iagemax+4+AGEMARGE; iage++) | for(iage=iagemin-AGEMARGE; iage <= iagemax+4+AGEMARGE; iage++) |
| Line 5088 void prevalence(double ***probs, double | Line 5157 void prevalence(double ***probs, double |
| if(posprop>=1.e-5){ | if(posprop>=1.e-5){ |
| probs[i][jk][j1]= prop[jk][i]/posprop; | probs[i][jk][j1]= prop[jk][i]/posprop; |
| } else{ | } else{ |
| if(first==1){ | if(!first){ |
| first=0; | first=1; |
| printf("Warning Observed prevalence doesn't sum to 1 for state %d: probs[%d][%d][%d]=%lf because of lack of cases\nSee others in log file...\n",jk,i,jk, j1,probs[i][jk][j1]); | printf("Warning Observed prevalence doesn't sum to 1 for state %d: probs[%d][%d][%d]=%lf because of lack of cases\nSee others in log file...\n",jk,i,jk, j1,probs[i][jk][j1]); |
| fprintf(ficlog,"Warning Observed prevalence doesn't sum to 1 for state %d: probs[%d][%d][%d]=%lf because of lack of cases\nSee others in log file...\n",jk,i,jk, j1,probs[i][jk][j1]); | |
| }else{ | }else{ |
| fprintf(ficlog,"Warning Observed prevalence doesn't sum to 1 for state %d: probs[%d][%d][%d]=%lf because of lack of cases\nSee others in log file...\n",jk,i,jk, j1,probs[i][jk][j1]); | fprintf(ficlog,"Warning Observed prevalence doesn't sum to 1 for state %d: probs[%d][%d][%d]=%lf because of lack of cases.\n",jk,i,jk, j1,probs[i][jk][j1]); |
| } | } |
| } | } |
| } | } |
| Line 5339 void concatwav(int wav[], int **dh, int | Line 5407 void concatwav(int wav[], int **dh, int |
| /* *cptcov=0; */ | /* *cptcov=0; */ |
| for (k=1; k <= maxncov; k++) ncodemax[k]=0; /* Horrible constant again replaced by NCOVMAX */ | for (k=1; k <= maxncov; k++) ncodemax[k]=0; /* Horrible constant again replaced by NCOVMAX */ |
| for (k=1; k <= maxncov; k++) | |
| for(j=1; j<=2; j++) | |
| nbcode[k][j]=0; /* Valgrind */ | |
| /* Loop on covariates without age and products and no quantitative variable */ | /* Loop on covariates without age and products and no quantitative variable */ |
| /* for (j=1; j<=(cptcovs); j++) { /\* From model V1 + V2*age+ V3 + V3*V4 keeps V1 + V3 = 2 only *\/ */ | |
| for (k=1; k<=cptcovt; k++) { /* From model V1 + V2*age + V3 + V3*V4 keeps V1 + V3 = 2 only */ | for (k=1; k<=cptcovt; k++) { /* From model V1 + V2*age + V3 + V3*V4 keeps V1 + V3 = 2 only */ |
| for (j=-1; (j < maxncov); j++) Ndum[j]=0; | for (j=-1; (j < maxncov); j++) Ndum[j]=0; |
| if(Dummy[k]==0 && Typevar[k] !=1){ /* Dummy covariate and not age product */ | if(Dummy[k]==0 && Typevar[k] !=1){ /* Dummy covariate and not age product */ |
| Line 5359 void concatwav(int wav[], int **dh, int | Line 5429 void concatwav(int wav[], int **dh, int |
| modmaxcovj=ij; | modmaxcovj=ij; |
| else if (ij < modmincovj) | else if (ij < modmincovj) |
| modmincovj=ij; | modmincovj=ij; |
| if ((ij < -1) && (ij > NCOVMAX)){ | if (ij <0 || ij >1 ){ |
| printf("Information, IMaCh doesn't treat covariate with missing values (-1), individual %d will be skipped.\n",i); | |
| fprintf(ficlog,"Information, currently IMaCh doesn't treat covariate with missing values (-1), individual %d will be skipped.\n",i); | |
| } | |
| if ((ij < -1) || (ij > NCOVMAX)){ | |
| printf( "Error: minimal is less than -1 or maximal is bigger than %d. Exiting. \n", NCOVMAX ); | printf( "Error: minimal is less than -1 or maximal is bigger than %d. Exiting. \n", NCOVMAX ); |
| exit(1); | exit(1); |
| }else | }else |
| Line 5405 void concatwav(int wav[], int **dh, int | Line 5479 void concatwav(int wav[], int **dh, int |
| /* nbcode[Tvar[j]][3]=2; */ | /* nbcode[Tvar[j]][3]=2; */ |
| /* To be continued (not working yet). */ | /* To be continued (not working yet). */ |
| ij=0; /* ij is similar to i but can jump over null modalities */ | ij=0; /* ij is similar to i but can jump over null modalities */ |
| for (i=modmincovj; i<=modmaxcovj; i++) { /* i= 1 to 2 for dichotomous, or from 1 to 3 or from -1 or 0 to 1 currently*/ | |
| /* for (i=modmincovj; i<=modmaxcovj; i++) { */ /* i= 1 to 2 for dichotomous, or from 1 to 3 or from -1 or 0 to 1 currently*/ | |
| /* Skipping the case of missing values by reducing nbcode to 0 and 1 and not -1, 0, 1 */ | |
| /* model=V1+V2+V3, if V2=-1, 0 or 1, then nbcode[2][1]=0 and nbcode[2][2]=1 instead of | |
| * nbcode[2][1]=-1, nbcode[2][2]=0 and nbcode[2][3]=1 */ | |
| /*, could be restored in the future */ | |
| for (i=0; i<=1; i++) { /* i= 1 to 2 for dichotomous, or from 1 to 3 or from -1 or 0 to 1 currently*/ | |
| if (Ndum[i] == 0) { /* If nobody responded to this modality k */ | if (Ndum[i] == 0) { /* If nobody responded to this modality k */ |
| break; | break; |
| } | } |
| ij++; | ij++; |
| nbcode[Tvar[k]][ij]=i; /* stores the original value of modality i in an array nbcode, ij modality from 1 to last non-nul modality. nbcode[1][1]=0 nbcode[1][2]=1*/ | nbcode[Tvar[k]][ij]=i; /* stores the original value of modality i in an array nbcode, ij modality from 1 to last non-nul modality. nbcode[1][1]=0 nbcode[1][2]=1 . Could be -1*/ |
| cptcode = ij; /* New max modality for covar j */ | cptcode = ij; /* New max modality for covar j */ |
| } /* end of loop on modality i=-1 to 1 or more */ | } /* end of loop on modality i=-1 to 1 or more */ |
| break; | break; |
| Line 5426 void concatwav(int wav[], int **dh, int | Line 5506 void concatwav(int wav[], int **dh, int |
| break; | break; |
| } /* end switch */ | } /* end switch */ |
| } /* end dummy test */ | } /* end dummy test */ |
| } /* end of loop on model-covariate k. nbcode[Tvark][1]=-1, nbcode[Tvark][1]=0 and nbcode[Tvark][2]=1 sets the value of covariate k*/ | |
| /* for (k=0; k<= cptcode; k++) { /\* k=-1 ? k=0 to 1 *\//\* Could be 1 to 4 *\//\* cptcode=modmaxcovj *\/ */ | |
| /* /\*recode from 0 *\/ */ | |
| /* k is a modality. If we have model=V1+V1*sex */ | |
| /* then: nbcode[1][1]=0 ; nbcode[1][2]=1; nbcode[2][1]=0 ; nbcode[2][2]=1; */ | |
| /* But if some modality were not used, it is recoded from 0 to a newer modmaxcovj=cptcode *\/ */ | |
| /* } */ | |
| /* /\* cptcode = ij; *\/ /\* New max modality for covar j *\/ */ | |
| /* if (ij > ncodemax[j]) { */ | |
| /* printf( " Error ij=%d > ncodemax[%d]=%d\n", ij, j, ncodemax[j]); */ | |
| /* fprintf(ficlog, " Error ij=%d > ncodemax[%d]=%d\n", ij, j, ncodemax[j]); */ | |
| /* break; */ | |
| /* } */ | |
| /* } /\* end of loop on modality k *\/ */ | |
| } /* end of loop on model-covariate j. nbcode[Tvarj][1]=0 and nbcode[Tvarj][2]=1 sets the value of covariate j*/ | |
| for (k=-1; k< maxncov; k++) Ndum[k]=0; | for (k=-1; k< maxncov; k++) Ndum[k]=0; |
| /* Look at fixed dummy (single or product) covariates to check empty modalities */ | /* Look at fixed dummy (single or product) covariates to check empty modalities */ |
| Line 5824 void concatwav(int wav[], int **dh, int | Line 5890 void concatwav(int wav[], int **dh, int |
| double **dnewm,**doldm; | double **dnewm,**doldm; |
| double **dnewmp,**doldmp; | double **dnewmp,**doldmp; |
| int i, j, nhstepm, hstepm, h, nstepm ; | int i, j, nhstepm, hstepm, h, nstepm ; |
| int first=0; | |
| int k; | int k; |
| double *xp; | double *xp; |
| double **gp, **gm; /**< for var eij */ | double **gp, **gm; /**< for var eij */ |
| Line 5948 void concatwav(int wav[], int **dh, int | Line 6015 void concatwav(int wav[], int **dh, int |
| } | } |
| /**< Computes the prevalence limit with parameter theta shifted of delta up to ftolpl precision and | /**< Computes the prevalence limit with parameter theta shifted of delta up to ftolpl precision and |
| * returns into prlim . | * returns into prlim . |
| */ | */ |
| prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij, nres); | prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij, nres); |
| /* If popbased = 1 we use crossection prevalences. Previous step is useless but prlim is created */ | /* If popbased = 1 we use crossection prevalences. Previous step is useless but prlim is created */ |
| Line 5964 void concatwav(int wav[], int **dh, int | Line 6031 void concatwav(int wav[], int **dh, int |
| /**< Computes the shifted transition matrix \f$ {}{h}_p^{ij}_x\f$ at horizon h. | /**< Computes the shifted transition matrix \f$ {}{h}_p^{ij}_x\f$ at horizon h. |
| */ | */ |
| hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij,nres); /* Returns p3mat[i][j][h] for h=0 to nhstepm */ | hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij,nres); /* Returns p3mat[i][j][h] for h=0 to nhstepm */ |
| /**< And for each alive state j, sums over i \f$ w^i_x {}{h}_p^{ij}_x\f$, which are the probability | /**< And for each alive state j, sums over i \f$ w^i_x {}{h}_p^{ij}x\f$, which are the probability |
| * at horizon h in state j including mortality. | * at horizon h in state j including mortality. |
| */ | */ |
| for(j=1; j<= nlstate; j++){ | for(j=1; j<= nlstate; j++){ |
| Line 5986 void concatwav(int wav[], int **dh, int | Line 6053 void concatwav(int wav[], int **dh, int |
| for(i=1; i<=npar; i++) /* Computes gradient x - delta */ | for(i=1; i<=npar; i++) /* Computes gradient x - delta */ |
| xp[i] = x[i] - (i==theta ?delti[theta]:0); | xp[i] = x[i] - (i==theta ?delti[theta]:0); |
| prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp, ij, nres); | prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp, ij, nres); |
| if (popbased==1) { | if (popbased==1) { |
| Line 6172 void concatwav(int wav[], int **dh, int | Line 6239 void concatwav(int wav[], int **dh, int |
| int theta; | int theta; |
| pstamp(ficresvpl); | pstamp(ficresvpl); |
| fprintf(ficresvpl,"# Standard deviation of period (stable) prevalences \n"); | fprintf(ficresvpl,"# Standard deviation of period (forward stable) prevalences \n"); |
| fprintf(ficresvpl,"# Age "); | fprintf(ficresvpl,"# Age "); |
| if(nresult >=1) | if(nresult >=1) |
| fprintf(ficresvpl," Result# "); | fprintf(ficresvpl," Result# "); |
| Line 6201 void concatwav(int wav[], int **dh, int | Line 6268 void concatwav(int wav[], int **dh, int |
| for(i=1; i<=npar; i++){ /* Computes gradient */ | for(i=1; i<=npar; i++){ /* Computes gradient */ |
| xp[i] = x[i] + (i==theta ?delti[theta]:0); | xp[i] = x[i] + (i==theta ?delti[theta]:0); |
| } | } |
| if((int)age==79 ||(int)age== 80 ||(int)age== 81 ) | /* if((int)age==79 ||(int)age== 80 ||(int)age== 81 ) */ |
| prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); | /* prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); */ |
| else | /* else */ |
| prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); | prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); |
| for(i=1;i<=nlstate;i++){ | for(i=1;i<=nlstate;i++){ |
| gp[i] = prlim[i][i]; | gp[i] = prlim[i][i]; |
| mgp[theta][i] = prlim[i][i]; | mgp[theta][i] = prlim[i][i]; |
| } | } |
| for(i=1; i<=npar; i++) /* Computes gradient */ | for(i=1; i<=npar; i++) /* Computes gradient */ |
| xp[i] = x[i] - (i==theta ?delti[theta]:0); | xp[i] = x[i] - (i==theta ?delti[theta]:0); |
| if((int)age==79 ||(int)age== 80 ||(int)age== 81 ) | /* if((int)age==79 ||(int)age== 80 ||(int)age== 81 ) */ |
| prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); | /* prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); */ |
| else | /* else */ |
| prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); | prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij,nres); |
| for(i=1;i<=nlstate;i++){ | for(i=1;i<=nlstate;i++){ |
| gm[i] = prlim[i][i]; | gm[i] = prlim[i][i]; |
| mgm[theta][i] = prlim[i][i]; | mgm[theta][i] = prlim[i][i]; |
| Line 6263 void concatwav(int wav[], int **dh, int | Line 6330 void concatwav(int wav[], int **dh, int |
| fprintf(ficresvpl,"%.0f ",age ); | fprintf(ficresvpl,"%.0f ",age ); |
| if(nresult >=1) | if(nresult >=1) |
| fprintf(ficresvpl,"%d ",nres ); | fprintf(ficresvpl,"%d ",nres ); |
| for(i=1; i<=nlstate;i++) | for(i=1; i<=nlstate;i++){ |
| fprintf(ficresvpl," %.5f (%.5f)",prlim[i][i],sqrt(varpl[i][(int)age])); | fprintf(ficresvpl," %.5f (%.5f)",prlim[i][i],sqrt(varpl[i][(int)age])); |
| /* for(j=1;j<=nlstate;j++) */ | |
| /* fprintf(ficresvpl," %d %.5f ",j,prlim[j][i]); */ | |
| } | |
| fprintf(ficresvpl,"\n"); | fprintf(ficresvpl,"\n"); |
| free_vector(gp,1,nlstate); | free_vector(gp,1,nlstate); |
| free_vector(gm,1,nlstate); | free_vector(gm,1,nlstate); |
| Line 6481 void varprob(char optionfilefiname[], do | Line 6551 void varprob(char optionfilefiname[], do |
| fprintf(fichtm,"\n<li><h4> Computing and drawing one step probabilities with their confidence intervals</h4></li>\n"); | fprintf(fichtm,"\n<li><h4> Computing and drawing one step probabilities with their confidence intervals</h4></li>\n"); |
| fprintf(fichtm,"\n"); | fprintf(fichtm,"\n"); |
| fprintf(fichtm,"\n<li><h4> <a href=\"%s\">Matrix of variance-covariance of one-step probabilities (drawings)</a></h4> this page is important in order to visualize confidence intervals and especially correlation between disability and recovery, or more generally, way in and way back. %s</li>\n",optionfilehtmcov,optionfilehtmcov); | fprintf(fichtm,"\n<li><h4> <a href=\"%s\">Matrix of variance-covariance of one-step probabilities (drawings)</a></h4> this page is important in order to visualize confidence intervals and especially correlation between disability and recovery, or more generally, way in and way back. File %s</li>\n",optionfilehtmcov,optionfilehtmcov); |
| fprintf(fichtmcov,"Current page is file <a href=\"%s\">%s</a><br>\n\n<h4>Matrix of variance-covariance of pairs of step probabilities</h4>\n",optionfilehtmcov, optionfilehtmcov); | fprintf(fichtmcov,"Current page is file <a href=\"%s\">%s</a><br>\n\n<h4>Matrix of variance-covariance of pairs of step probabilities</h4>\n",optionfilehtmcov, optionfilehtmcov); |
| fprintf(fichtmcov,"\nEllipsoids of confidence centered on point (p<inf>ij</inf>, p<inf>kl</inf>) are estimated \ | fprintf(fichtmcov,"\nEllipsoids of confidence centered on point (p<inf>ij</inf>, p<inf>kl</inf>) are estimated \ |
| and drawn. It helps understanding how is the covariance between two incidences.\ | and drawn. It helps understanding how is the covariance between two incidences.\ |
| Line 6774 void printinghtml(char fileresu[], char | Line 6844 void printinghtml(char fileresu[], char |
| - Estimated back transition probabilities over %d (stepm) months: <a href=\"%s\">%s</a><br>\n ", | - Estimated back transition probabilities over %d (stepm) months: <a href=\"%s\">%s</a><br>\n ", |
| stepm,subdirf2(fileresu,"PIJB_"),subdirf2(fileresu,"PIJB_")); | stepm,subdirf2(fileresu,"PIJB_"),subdirf2(fileresu,"PIJB_")); |
| fprintf(fichtm,"\ | fprintf(fichtm,"\ |
| - Period (stable) prevalence in each health state: <a href=\"%s\">%s</a> <br>\n", | - Period (forward) prevalence in each health state: <a href=\"%s\">%s</a> <br>\n", |
| subdirf2(fileresu,"PL_"),subdirf2(fileresu,"PL_")); | subdirf2(fileresu,"PL_"),subdirf2(fileresu,"PL_")); |
| fprintf(fichtm,"\ | fprintf(fichtm,"\ |
| - Period (stable) back prevalence in each health state: <a href=\"%s\">%s</a> <br>\n", | - Backward prevalence in each health state: <a href=\"%s\">%s</a> <br>\n", |
| subdirf2(fileresu,"PLB_"),subdirf2(fileresu,"PLB_")); | subdirf2(fileresu,"PLB_"),subdirf2(fileresu,"PLB_")); |
| fprintf(fichtm,"\ | fprintf(fichtm,"\ |
| - (a) Life expectancies by health status at initial age, e<sub>i.</sub> (b) health expectancies by health status at initial age, e<sub>ij</sub> . If one or more covariates are included, specific tables for each value of the covariate are output in sequences within the same file (estepm=%2d months): \ | - (a) Life expectancies by health status at initial age, e<sub>i.</sub> (b) health expectancies by health status at initial age, e<sub>ij</sub> . If one or more covariates are included, specific tables for each value of the covariate are output in sequences within the same file (estepm=%2d months): \ |
| Line 6883 divided by h: <sub>h</sub>P<sub>ij</sub> | Line 6953 divided by h: <sub>h</sub>P<sub>ij</sub> |
| <img src=\"%s_%d-3-%d.svg\">",stepm,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres); | <img src=\"%s_%d-3-%d.svg\">",stepm,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres,subdirf2(optionfilefiname,"PE_"),k1,nres); |
| /* Survival functions (period) in state j */ | /* Survival functions (period) in state j */ |
| for(cpt=1; cpt<=nlstate;cpt++){ | for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Survival functions in state %d. Or probability to survive in state %d being in state (1 to %d) at different ages. <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ | fprintf(fichtm,"<br>\n- Survival functions in state %d. And probability to be observed in state %d being in state (1 to %d) at different ages. <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ |
| <img src=\"%s_%d-%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres); | <img src=\"%s_%d-%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJ_"),cpt,k1,nres); |
| } | } |
| /* State specific survival functions (period) */ | /* State specific survival functions (period) */ |
| for(cpt=1; cpt<=nlstate;cpt++){ | for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Survival functions from state %d in each live state and total.\ | fprintf(fichtm,"<br>\n- Survival functions in state %d and in any other live state (total).\ |
| Or probability to survive in various states (1 to %d) being in state %d at different ages. \ | And probability to be observed in various states (up to %d) being in state %d at different ages. \ |
| <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> <img src=\"%s_%d-%d-%d.svg\">", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres); | <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> <img src=\"%s_%d-%d-%d.svg\">", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres); |
| } | } |
| /* Period (stable) prevalence in each health state */ | /* Period (forward stable) prevalence in each health state */ |
| for(cpt=1; cpt<=nlstate;cpt++){ | for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Convergence to period (stable) prevalence in state %d. Or probability for a person being in state (1 to %d) at different ages, to be in state %d some years after. <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ | fprintf(fichtm,"<br>\n- Convergence to period (stable) prevalence in state %d. Or probability for a person being in state (1 to %d) at different ages, to be in state %d some years after. <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ |
| <img src=\"%s_%d-%d-%d.svg\">", cpt, nlstate, cpt, subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres); | <img src=\"%s_%d-%d-%d.svg\">", cpt, nlstate, cpt, subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres,subdirf2(optionfilefiname,"P_"),cpt,k1,nres); |
| } | } |
| if(backcast==1){ | if(backcast==1){ |
| /* Period (stable) back prevalence in each health state */ | /* Backward prevalence in each health state */ |
| for(cpt=1; cpt<=nlstate;cpt++){ | for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Convergence to mixed (stable) back prevalence in state %d. Or probability for a person to be in state %d at a younger age, knowing that she/he was in state (1 to %d) at different older ages. <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ | fprintf(fichtm,"<br>\n- Convergence to mixed (stable) back prevalence in state %d. Or probability for a person to be in state %d at a younger age, knowing that she/he was in state (1 to %d) at different older ages. <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ |
| <img src=\"%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); | <img src=\"%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){ | if(prevfcast==1){ |
| /* Projection of prevalence up to period (stable) prevalence in each health state */ | /* Projection of prevalence up to period (forward stable) prevalence in each health state */ |
| for(cpt=1; cpt<=nlstate;cpt++){ | for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f and mobil_average=%d), from year %.1f up to year %.1f tending to period (stable) prevalence in state %d. Or probability to be in state %d being in an observed weighted state (from 1 to %d). <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ | fprintf(fichtm,"<br>\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f and mobil_average=%d), from year %.1f up to year %.1f tending to period (stable) forward prevalence in state %d. Or probability to be in state %d being in an observed weighted state (from 1 to %d). <a href=\"%s_%d-%d-%d.svg\">%s_%d-%d-%d.svg</a><br> \ |
| <img src=\"%s_%d-%d-%d.svg\">", dateprev1, dateprev2, mobilavproj, dateproj1, dateproj2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres); | <img src=\"%s_%d-%d-%d.svg\">", dateprev1, dateprev2, mobilavproj, dateproj1, dateproj2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres,subdirf2(optionfilefiname,"PROJ_"),cpt,k1,nres); |
| } | } |
| } | } |
| Line 6960 See page 'Matrix of variance-covariance | Line 7030 See page 'Matrix of variance-covariance |
| <a href=\"%s\">%s</a> <br>\n</li>", | <a href=\"%s\">%s</a> <br>\n</li>", |
| estepm,subdirf2(fileresu,"STDE_"),subdirf2(fileresu,"STDE_")); | estepm,subdirf2(fileresu,"STDE_"),subdirf2(fileresu,"STDE_")); |
| fprintf(fichtm,"\ | fprintf(fichtm,"\ |
| - Variances and covariances of health expectancies by age. Status (i) based health expectancies (in state j), e<sup>ij</sup> are weighted by the period prevalences in each state i (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a><br>\n", | - Variances and covariances of health expectancies by age. Status (i) based health expectancies (in state j), e<sup>ij</sup> are weighted by the forward (period) prevalences in each state i (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a><br>\n", |
| estepm, subdirf2(fileresu,"V_"),subdirf2(fileresu,"V_")); | estepm, subdirf2(fileresu,"V_"),subdirf2(fileresu,"V_")); |
| fprintf(fichtm,"\ | fprintf(fichtm,"\ |
| - Total life expectancy and total health expectancies to be spent in each health state e<sup>.j</sup> with their standard errors (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a> <br>\n", | - Total life expectancy and total health expectancies to be spent in each health state e<sup>.j</sup> with their standard errors (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a> <br>\n", |
| estepm, subdirf2(fileresu,"T_"),subdirf2(fileresu,"T_")); | estepm, subdirf2(fileresu,"T_"),subdirf2(fileresu,"T_")); |
| fprintf(fichtm,"\ | fprintf(fichtm,"\ |
| - Standard deviation of period (stable) prevalences: <a href=\"%s\">%s</a> <br>\n",\ | - Standard deviation of forward (period) prevalences: <a href=\"%s\">%s</a> <br>\n",\ |
| subdirf2(fileresu,"VPL_"),subdirf2(fileresu,"VPL_")); | subdirf2(fileresu,"VPL_"),subdirf2(fileresu,"VPL_")); |
| /* if(popforecast==1) fprintf(fichtm,"\n */ | /* if(popforecast==1) fprintf(fichtm,"\n */ |
| Line 7102 void printinggnuplot(char fileresu[], ch | Line 7172 void printinggnuplot(char fileresu[], ch |
| continue; | continue; |
| /* We are interested in selected combination by the resultline */ | /* We are interested in selected combination by the resultline */ |
| /* printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt); */ | /* printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt); */ |
| fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt); | fprintf(ficgp,"\n# 1st: Forward (stable period) prevalence with CI: 'VPL_' files and live state =%d ", cpt); |
| strcpy(gplotlabel,"("); | strcpy(gplotlabel,"("); |
| for (k=1; k<=cptcoveff; k++){ /* For each covariate k get corresponding value lv for combination k1 */ | for (k=1; k<=cptcoveff; k++){ /* For each covariate k get corresponding value lv for combination k1 */ |
| lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */ | lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */ |
| Line 7140 void printinggnuplot(char fileresu[], ch | Line 7210 void printinggnuplot(char fileresu[], ch |
| if (i==cpt) fprintf(ficgp," %%lf (%%lf)"); | if (i==cpt) fprintf(ficgp," %%lf (%%lf)"); |
| else fprintf(ficgp," %%*lf (%%*lf)"); | else fprintf(ficgp," %%*lf (%%*lf)"); |
| } | } |
| fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2==%d ? $3+1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),nres-1,nres-1,nres); | fprintf(ficgp,"\" t\"Forward prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2==%d ? $3+1.96*$4 : 1/0) \"%%lf %%lf",subdirf2(fileresu,"VPL_"),nres-1,nres-1,nres); |
| for (i=1; i<= nlstate ; i ++) { | for (i=1; i<= nlstate ; i ++) { |
| if (i==cpt) fprintf(ficgp," %%lf (%%lf)"); | if (i==cpt) fprintf(ficgp," %%lf (%%lf)"); |
| else fprintf(ficgp," %%*lf (%%*lf)"); | else fprintf(ficgp," %%*lf (%%*lf)"); |
| Line 7469 set ter svg size 640, 480\nunset log y\n | Line 7539 set ter svg size 640, 480\nunset log y\n |
| continue; | continue; |
| for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state of arrival */ | for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state of arrival */ |
| strcpy(gplotlabel,"("); | strcpy(gplotlabel,"("); |
| fprintf(ficgp,"\n#\n#\n#CV preval stable (period): 'pij' files, covariatecombination#=%d state=%d",k1, cpt); | fprintf(ficgp,"\n#\n#\n#CV preval stable (forward): 'pij' files, covariatecombination#=%d state=%d",k1, cpt); |
| for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */ | 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 */ | 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 */ | /* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */ |
| Line 7513 set ter svg size 640, 480\nunset log y\n | Line 7583 set ter svg size 640, 480\nunset log y\n |
| /* 7eme */ | /* 7eme */ |
| if(backcast == 1){ | if(backcast == 1){ |
| /* CV back preval stable (period) for each covariate */ | /* CV backward prevalence for each covariate */ |
| for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ | for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ |
| for(nres=1; nres <= nresult; nres++){ /* For each resultline */ | for(nres=1; nres <= nresult; nres++){ /* For each resultline */ |
| if(m != 1 && TKresult[nres]!= k1) | if(m != 1 && TKresult[nres]!= k1) |
| continue; | continue; |
| for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life origin state */ | for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life origin state */ |
| strcpy(gplotlabel,"("); | strcpy(gplotlabel,"("); |
| fprintf(ficgp,"\n#\n#\n#CV Back preval stable (period): 'pijb' files, covariatecombination#=%d state=%d",k1, cpt); | fprintf(ficgp,"\n#\n#\n#CV Backward stable prevalence: 'pijb' files, covariatecombination#=%d state=%d",k1, cpt); |
| for (k=1; k<=cptcoveff; k++){ /* For each covariate and each value */ | 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 */ | 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 */ | /* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */ |
| Line 7568 set ter svg size 640, 480\nunset log y\n | Line 7638 set ter svg size 640, 480\nunset log y\n |
| /* 8eme */ | /* 8eme */ |
| if(prevfcast==1){ | if(prevfcast==1){ |
| /* Projection from cross-sectional to stable (period) for each covariate */ | /* Projection from cross-sectional to forward stable (period) prevalence for each covariate */ |
| for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ | for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */ |
| for(nres=1; nres <= nresult; nres++){ /* For each resultline */ | for(nres=1; nres <= nresult; nres++){ /* For each resultline */ |
| Line 7576 set ter svg size 640, 480\nunset log y\n | Line 7646 set ter svg size 640, 480\nunset log y\n |
| continue; | continue; |
| for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */ | for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */ |
| strcpy(gplotlabel,"("); | strcpy(gplotlabel,"("); |
| fprintf(ficgp,"\n#\n#\n#Projection of prevalence to stable (period): 'PROJ_' files, covariatecombination#=%d state=%d",k1, cpt); | fprintf(ficgp,"\n#\n#\n#Projection of prevalence to forward stable prevalence (period): 'PROJ_' files, covariatecombination#=%d state=%d",k1, cpt); |
| for (k=1; k<=cptcoveff; k++){ /* For each correspondig covariate value */ | for (k=1; k<=cptcoveff; k++){ /* For each correspondig covariate value */ |
| lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate value corresponding to k1 combination and kth covariate */ | lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate value corresponding to k1 combination and kth covariate */ |
| /* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */ | /* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */ |
| Line 8010 set ter svg size 640, 480\nunset log y\n | Line 8080 set ter svg size 640, 480\nunset log y\n |
| int modcovmax =1; | int modcovmax =1; |
| int mobilavrange, mob; | int mobilavrange, mob; |
| int iage=0; | int iage=0; |
| int firstA1=0, firstA2=0; | |
| double sum=0., sumr=0.; | double sum=0., sumr=0.; |
| double age; | double age; |
| Line 8107 set ter svg size 640, 480\nunset log y\n | Line 8178 set ter svg size 640, 480\nunset log y\n |
| } /* age */ | } /* age */ |
| /* Thus we have agemingood and agemaxgood as well as goodr for raw (preobs) */ | /* Thus we have agemingood and agemaxgood as well as goodr for raw (preobs) */ |
| /* but they will change */ | /* but they will change */ |
| firstA1=0;firstA2=0; | |
| for (age=fage-(mob-1)/2; age>=bage; age--){/* From oldest to youngest, filling up to the youngest */ | for (age=fage-(mob-1)/2; age>=bage; age--){/* From oldest to youngest, filling up to the youngest */ |
| sumnewm[cptcod]=0.; | sumnewm[cptcod]=0.; |
| sumnewmr[cptcod]=0.; | sumnewmr[cptcod]=0.; |
| Line 8139 set ter svg size 640, 480\nunset log y\n | Line 8211 set ter svg size 640, 480\nunset log y\n |
| sumr+=probs[(int)age][i][cptcod]; | sumr+=probs[(int)age][i][cptcod]; |
| } | } |
| if(fabs(sum - 1.) > 1.e-3) { /* bad */ | if(fabs(sum - 1.) > 1.e-3) { /* bad */ |
| printf("Moving average A1: For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one (%f) at any descending age! age=%d, could you increase bage=%d\n",cptcod,sumr, (int)age, (int)bage); | if(!firstA1){ |
| firstA1=1; | |
| printf("Moving average A1: For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one (%f) at any descending age! age=%d, could you increase bage=%d. Others in log file...\n",cptcod,sumr, (int)age, (int)bage); | |
| } | |
| fprintf(ficlog,"Moving average A1: For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one (%f) at any descending age! age=%d, could you increase bage=%d\n",cptcod,sumr, (int)age, (int)bage); | |
| } /* end bad */ | } /* end bad */ |
| /* else{ /\* We found some ages summing to one, we will smooth the oldest *\/ */ | /* else{ /\* We found some ages summing to one, we will smooth the oldest *\/ */ |
| if(fabs(sumr - 1.) > 1.e-3) { /* bad */ | if(fabs(sumr - 1.) > 1.e-3) { /* bad */ |
| printf("Moving average A2: For this combination of covariate cptcod=%d, the raw prevalence doesn't sums to one (%f) even with smoothed values at young ages! age=%d, could you increase bage=%d\n",cptcod,sumr, (int)age, (int)bage); | if(!firstA2){ |
| firstA2=1; | |
| printf("Moving average A2: For this combination of covariate cptcod=%d, the raw prevalence doesn't sums to one (%f) even with smoothed values at young ages! age=%d, could you increase bage=%d. Others in log file...\n",cptcod,sumr, (int)age, (int)bage); | |
| } | |
| fprintf(ficlog,"Moving average A2: For this combination of covariate cptcod=%d, the raw prevalence doesn't sums to one (%f) even with smoothed values at young ages! age=%d, could you increase bage=%d\n",cptcod,sumr, (int)age, (int)bage); | |
| } /* end bad */ | } /* end bad */ |
| }/* age */ | }/* age */ |
| Line 8518 set ter svg size 640, 480\nunset log y\n | Line 8598 set ter svg size 640, 480\nunset log y\n |
| /* Variance of prevalence limit: varprlim */ | /* Variance of prevalence limit: varprlim */ |
| void varprlim(char fileresu[], int nresult, double ***prevacurrent, int mobilavproj, double bage, double fage, double **prlim, int *ncvyearp, double ftolpl, double p[], double **matcov, double *delti, int stepm, int cptcoveff){ | void varprlim(char fileresu[], int nresult, double ***prevacurrent, int mobilavproj, double bage, double fage, double **prlim, int *ncvyearp, double ftolpl, double p[], double **matcov, double *delti, int stepm, int cptcoveff){ |
| /*------- Variance of period (stable) prevalence------*/ | /*------- Variance of forward period (stable) prevalence------*/ |
| char fileresvpl[FILENAMELENGTH]; | char fileresvpl[FILENAMELENGTH]; |
| FILE *ficresvpl; | FILE *ficresvpl; |
| Line 8529 set ter svg size 640, 480\nunset log y\n | Line 8609 set ter svg size 640, 480\nunset log y\n |
| strcpy(fileresvpl,"VPL_"); | strcpy(fileresvpl,"VPL_"); |
| strcat(fileresvpl,fileresu); | strcat(fileresvpl,fileresu); |
| if((ficresvpl=fopen(fileresvpl,"w"))==NULL) { | if((ficresvpl=fopen(fileresvpl,"w"))==NULL) { |
| printf("Problem with variance of period (stable) prevalence resultfile: %s\n", fileresvpl); | printf("Problem with variance of forward period (stable) prevalence resultfile: %s\n", fileresvpl); |
| exit(0); | exit(0); |
| } | } |
| printf("Computing Variance-covariance of period (stable) prevalence: file '%s' ...", fileresvpl);fflush(stdout); | printf("Computing Variance-covariance of forward period (stable) prevalence: file '%s' ...", fileresvpl);fflush(stdout); |
| fprintf(ficlog, "Computing Variance-covariance of period (stable) prevalence: file '%s' ...", fileresvpl);fflush(ficlog); | fprintf(ficlog, "Computing Variance-covariance of forward period (stable) prevalence: file '%s' ...", fileresvpl);fflush(ficlog); |
| /*for(cptcov=1,k=0;cptcov<=i1;cptcov++){ | /*for(cptcov=1,k=0;cptcov<=i1;cptcov++){ |
| for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*/ | for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*/ |
| Line 8570 set ter svg size 640, 480\nunset log y\n | Line 8650 set ter svg size 640, 480\nunset log y\n |
| } | } |
| fclose(ficresvpl); | fclose(ficresvpl); |
| printf("done variance-covariance of period prevalence\n");fflush(stdout); | printf("done variance-covariance of forward period prevalence\n");fflush(stdout); |
| fprintf(ficlog,"done variance-covariance of period prevalence\n");fflush(ficlog); | fprintf(ficlog,"done variance-covariance of forward period prevalence\n");fflush(ficlog); |
| } | } |
| /* Variance of back prevalence: varbprlim */ | /* Variance of back prevalence: varbprlim */ |
| Line 9714 Dummy[k] 0=dummy (0 1), 1 quantitative ( | Line 9794 Dummy[k] 0=dummy (0 1), 1 quantitative ( |
| Typevar: 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product \n\ | Typevar: 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for product \n\ |
| Fixed[k] 0=fixed (product or simple), 1 varying, 2 fixed with age product, 3 varying with age product \n\ | Fixed[k] 0=fixed (product or simple), 1 varying, 2 fixed with age product, 3 varying with age product \n\ |
| Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product\n",model); | Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product\n",model); |
| for(k=1;k<=cptcovt; k++){ Fixed[k]=0; Dummy[k]=0;} | for(k=-1;k<=cptcovt; k++){ Fixed[k]=0; Dummy[k]=0;} |
| for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0;k<=cptcovt; k++){ /* or cptocvt */ | for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0;k<=cptcovt; k++){ /* or cptocvt */ |
| if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */ | if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */ |
| Fixed[k]= 0; | Fixed[k]= 0; |
| Line 9964 Dummy[k] 0=dummy (0 1), 1 quantitative ( | Line 10044 Dummy[k] 0=dummy (0 1), 1 quantitative ( |
| /* Searching for doublons in the model */ | /* Searching for doublons in the model */ |
| for(k1=1; k1<= cptcovt;k1++){ | for(k1=1; k1<= cptcovt;k1++){ |
| for(k2=1; k2 <k1;k2++){ | for(k2=1; k2 <k1;k2++){ |
| if((Typevar[k1]==Typevar[k2]) && (Fixed[Tvar[k1]]==Fixed[Tvar[k2]]) && (Dummy[Tvar[k1]]==Dummy[Tvar[k2]] )){ | /* if((Typevar[k1]==Typevar[k2]) && (Fixed[Tvar[k1]]==Fixed[Tvar[k2]]) && (Dummy[Tvar[k1]]==Dummy[Tvar[k2]] )){ */ |
| if((Typevar[k1]==Typevar[k2]) && (Fixed[k1]==Fixed[k2]) && (Dummy[k1]==Dummy[k2] )){ | |
| if((Typevar[k1] == 0 || Typevar[k1] == 1)){ /* Simple or age product */ | if((Typevar[k1] == 0 || Typevar[k1] == 1)){ /* Simple or age product */ |
| if(Tvar[k1]==Tvar[k2]){ | if(Tvar[k1]==Tvar[k2]){ |
| printf("Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); | printf("Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[k1],Dummy[k1]); |
| fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[Tvar[k1]],Dummy[Tvar[k1]]); fflush(ficlog); | fprintf(ficlog,"Error duplication in the model=%s at positions (+) %d and %d, Tvar[%d]=V%d, Tvar[%d]=V%d, Typevar=%d, Fixed=%d, Dummy=%d\n", model, k1,k2, k1, Tvar[k1], k2, Tvar[k2],Typevar[k1],Fixed[k1],Dummy[k1]); fflush(ficlog); |
| return(1); | return(1); |
| } | } |
| }else if (Typevar[k1] ==2){ | }else if (Typevar[k1] ==2){ |
| Line 10140 BOOL IsWow64() | Line 10221 BOOL IsWow64() |
| #endif | #endif |
| void syscompilerinfo(int logged) | void syscompilerinfo(int logged) |
| { | { |
| /* #include "syscompilerinfo.h"*/ | #include <stdint.h> |
| /* #include "syscompilerinfo.h"*/ | |
| /* command line Intel compiler 32bit windows, XP compatible:*/ | /* command line Intel compiler 32bit windows, XP compatible:*/ |
| /* /GS /W3 /Gy | /* /GS /W3 /Gy |
| /Zc:wchar_t /Zi /O2 /Fd"Release\vc120.pdb" /D "WIN32" /D "NDEBUG" /D | /Zc:wchar_t /Zi /O2 /Fd"Release\vc120.pdb" /D "WIN32" /D "NDEBUG" /D |
| Line 10176 void syscompilerinfo(int logged) | Line 10259 void syscompilerinfo(int logged) |
| /ManifestFile:"Release\IMaCh.exe.intermediate.manifest" /OPT:ICF | /ManifestFile:"Release\IMaCh.exe.intermediate.manifest" /OPT:ICF |
| /NOLOGO /TLBID:1 | /NOLOGO /TLBID:1 |
| */ | */ |
| #if defined __INTEL_COMPILER | #if defined __INTEL_COMPILER |
| #if defined(__GNUC__) | #if defined(__GNUC__) |
| struct utsname sysInfo; /* For Intel on Linux and OS/X */ | struct utsname sysInfo; /* For Intel on Linux and OS/X */ |
| Line 10192 void syscompilerinfo(int logged) | Line 10277 void syscompilerinfo(int logged) |
| } | } |
| #endif | #endif |
| #include <stdint.h> | |
| printf("Compiled with:");if(logged)fprintf(ficlog,"Compiled with:"); | printf("Compiled with:");if(logged)fprintf(ficlog,"Compiled with:"); |
| #if defined(__clang__) | #if defined(__clang__) |
| printf(" Clang/LLVM");if(logged)fprintf(ficlog," Clang/LLVM"); /* Clang/LLVM. ---------------------------------------------- */ | printf(" Clang/LLVM");if(logged)fprintf(ficlog," Clang/LLVM"); /* Clang/LLVM. ---------------------------------------------- */ |
| Line 10279 void syscompilerinfo(int logged) | Line 10362 void syscompilerinfo(int logged) |
| #endif | #endif |
| #endif | #endif |
| // void main() | // void main () |
| // { | // { |
| #if defined(_MSC_VER) | #if defined(_MSC_VER) |
| if (IsWow64()){ | if (IsWow64()){ |
| Line 10300 void syscompilerinfo(int logged) | Line 10383 void syscompilerinfo(int logged) |
| } | } |
| int prevalence_limit(double *p, double **prlim, double ageminpar, double agemaxpar, double ftolpl, int *ncvyearp){ | int prevalence_limit(double *p, double **prlim, double ageminpar, double agemaxpar, double ftolpl, int *ncvyearp){ |
| /*--------------- Prevalence limit (period or stable prevalence) --------------*/ | /*--------------- Prevalence limit (forward period or forward stable prevalence) --------------*/ |
| int i, j, k, i1, k4=0, nres=0 ; | int i, j, k, i1, k4=0, nres=0 ; |
| /* double ftolpl = 1.e-10; */ | /* double ftolpl = 1.e-10; */ |
| double age, agebase, agelim; | double age, agebase, agelim; |
| Line 10309 int prevalence_limit(double *p, double * | Line 10392 int prevalence_limit(double *p, double * |
| strcpy(filerespl,"PL_"); | strcpy(filerespl,"PL_"); |
| strcat(filerespl,fileresu); | strcat(filerespl,fileresu); |
| if((ficrespl=fopen(filerespl,"w"))==NULL) { | if((ficrespl=fopen(filerespl,"w"))==NULL) { |
| printf("Problem with period (stable) prevalence resultfile: %s\n", filerespl);return 1; | printf("Problem with forward period (stable) prevalence resultfile: %s\n", filerespl);return 1; |
| fprintf(ficlog,"Problem with period (stable) prevalence resultfile: %s\n", filerespl);return 1; | fprintf(ficlog,"Problem with forward period (stable) prevalence resultfile: %s\n", filerespl);return 1; |
| } | } |
| printf("\nComputing period (stable) prevalence: result on file '%s' \n", filerespl); | printf("\nComputing forward period (stable) prevalence: result on file '%s' \n", filerespl); |
| fprintf(ficlog,"\nComputing period (stable) prevalence: result on file '%s' \n", filerespl); | fprintf(ficlog,"\nComputing forward period (stable) prevalence: result on file '%s' \n", filerespl); |
| pstamp(ficrespl); | pstamp(ficrespl); |
| fprintf(ficrespl,"# Period (stable) prevalence. Precision given by ftolpl=%g \n", ftolpl); | fprintf(ficrespl,"# Forward period (stable) prevalence. Precision given by ftolpl=%g \n", ftolpl); |
| fprintf(ficrespl,"#Age "); | fprintf(ficrespl,"#Age "); |
| for(i=1; i<=nlstate;i++) fprintf(ficrespl,"%d-%d ",i,i); | for(i=1; i<=nlstate;i++) fprintf(ficrespl,"%d-%d ",i,i); |
| fprintf(ficrespl,"\n"); | fprintf(ficrespl,"\n"); |
| Line 10390 int prevalence_limit(double *p, double * | Line 10473 int prevalence_limit(double *p, double * |
| } | } |
| int back_prevalence_limit(double *p, double **bprlim, double ageminpar, double agemaxpar, double ftolpl, int *ncvyearp, double dateprev1,double dateprev2, int firstpass, int lastpass, int mobilavproj){ | int back_prevalence_limit(double *p, double **bprlim, double ageminpar, double agemaxpar, double ftolpl, int *ncvyearp, double dateprev1,double dateprev2, int firstpass, int lastpass, int mobilavproj){ |
| /*--------------- Back Prevalence limit (period or stable prevalence) --------------*/ | /*--------------- Back Prevalence limit (backward stable prevalence) --------------*/ |
| /* Computes the back prevalence limit for any combination of covariate values | /* Computes the back prevalence limit for any combination of covariate values |
| * at any age between ageminpar and agemaxpar | * at any age between ageminpar and agemaxpar |
| Line 10405 int back_prevalence_limit(double *p, dou | Line 10488 int back_prevalence_limit(double *p, dou |
| strcpy(fileresplb,"PLB_"); | strcpy(fileresplb,"PLB_"); |
| strcat(fileresplb,fileresu); | strcat(fileresplb,fileresu); |
| if((ficresplb=fopen(fileresplb,"w"))==NULL) { | if((ficresplb=fopen(fileresplb,"w"))==NULL) { |
| printf("Problem with period (stable) back prevalence resultfile: %s\n", fileresplb);return 1; | printf("Problem with backward prevalence resultfile: %s\n", fileresplb);return 1; |
| fprintf(ficlog,"Problem with period (stable) back prevalence resultfile: %s\n", fileresplb);return 1; | fprintf(ficlog,"Problem with backward prevalence resultfile: %s\n", fileresplb);return 1; |
| } | } |
| printf("Computing period (stable) back prevalence: result on file '%s' \n", fileresplb); | printf("Computing backward prevalence: result on file '%s' \n", fileresplb); |
| fprintf(ficlog,"Computing period (stable) back prevalence: result on file '%s' \n", fileresplb); | fprintf(ficlog,"Computing backward prevalence: result on file '%s' \n", fileresplb); |
| pstamp(ficresplb); | pstamp(ficresplb); |
| fprintf(ficresplb,"# Period (stable) back prevalence. Precision given by ftolpl=%g \n", ftolpl); | fprintf(ficresplb,"# Backward prevalence. Precision given by ftolpl=%g \n", ftolpl); |
| fprintf(ficresplb,"#Age "); | fprintf(ficresplb,"#Age "); |
| for(i=1; i<=nlstate;i++) fprintf(ficresplb,"%d-%d ",i,i); | for(i=1; i<=nlstate;i++) fprintf(ficresplb,"%d-%d ",i,i); |
| fprintf(ficresplb,"\n"); | fprintf(ficresplb,"\n"); |
| Line 10600 int hPijx(double *p, int bage, int fage) | Line 10683 int hPijx(double *p, int bage, int fage) |
| /*if (stepm<=24) stepsize=2;*/ | /*if (stepm<=24) stepsize=2;*/ |
| /* agelim=AGESUP; */ | /* agelim=AGESUP; */ |
| ageminl=30; | ageminl=AGEINF; /* was 30 */ |
| hstepm=stepsize*YEARM; /* Every year of age */ | hstepm=stepsize*YEARM; /* Every year of age */ |
| hstepm=hstepm/stepm; /* Typically 2 years, = 2/6 months = 4 */ | hstepm=hstepm/stepm; /* Typically 2 years, = 2/6 months = 4 */ |
| Line 10679 int main(int argc, char *argv[]) | Line 10762 int main(int argc, char *argv[]) |
| double ssval; | double ssval; |
| #endif | #endif |
| int movingaverage(double ***probs, double bage,double fage, double ***mobaverage, int mobilav); | int movingaverage(double ***probs, double bage,double fage, double ***mobaverage, int mobilav); |
| int i,j, k, n=MAXN,iter=0,m,size=100, cptcod; | int i,j, k, iter=0,m,size=100, cptcod; /* Suppressing because nobs */ |
| /* int i,j, k, n=MAXN,iter=0,m,size=100, cptcod; */ | |
| int ncvyear=0; /* Number of years needed for the period prevalence to converge */ | int ncvyear=0; /* Number of years needed for the period prevalence to converge */ |
| int jj, ll, li, lj, lk; | int jj, ll, li, lj, lk; |
| int numlinepar=0; /* Current linenumber of parameter file */ | int numlinepar=0; /* Current linenumber of parameter file */ |
| Line 10714 int main(int argc, char *argv[]) | Line 10798 int main(int argc, char *argv[]) |
| char pathr[MAXLINE], pathimach[MAXLINE]; | char pathr[MAXLINE], pathimach[MAXLINE]; |
| char *tok, *val; /* pathtot */ | char *tok, *val; /* pathtot */ |
| int firstobs=1, lastobs=10; | int firstobs=1, lastobs=10; /* nobs = lastobs-firstobs declared globally ;*/ |
| int c, h , cpt, c2; | int c, h , cpt, c2; |
| int jl=0; | int jl=0; |
| int i1, j1, jk, stepsize=0; | int i1, j1, jk, stepsize=0; |
| Line 10722 int main(int argc, char *argv[]) | Line 10806 int main(int argc, char *argv[]) |
| int *tab; | int *tab; |
| int mobilavproj=0 , prevfcast=0 ; /* moving average of prev, If prevfcast=1 prevalence projection */ | int mobilavproj=0 , prevfcast=0 ; /* moving average of prev, If prevfcast=1 prevalence projection */ |
| int backcast=0; | int backcast=0; /* defined as global for mlikeli and mle*/ |
| int mobilav=0,popforecast=0; | int mobilav=0,popforecast=0; |
| int hstepm=0, nhstepm=0; | int hstepm=0, nhstepm=0; |
| int agemortsup; | int agemortsup; |
| Line 11003 int main(int argc, char *argv[]) | Line 11087 int main(int argc, char *argv[]) |
| fprintf(ficlog,"Not 11 parameters, for example:ftol=1.e-8 stepm=12 ncovcol=2 nqv=1 ntv=2 nqtv=1 nlstate=2 ndeath=1 maxwav=3 mle=1 weight=1\n"); | fprintf(ficlog,"Not 11 parameters, for example:ftol=1.e-8 stepm=12 ncovcol=2 nqv=1 ntv=2 nqtv=1 nlstate=2 ndeath=1 maxwav=3 mle=1 weight=1\n"); |
| fprintf(ficlog,"but line=%s\n",line); | fprintf(ficlog,"but line=%s\n",line); |
| } | } |
| printf("ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); | if( lastpass > maxwav){ |
| printf("Error (lastpass = %d) > (maxwav = %d)\n",lastpass, maxwav); | |
| fprintf(ficlog,"Error (lastpass = %d) > (maxwav = %d)\n",lastpass, maxwav); | |
| fflush(ficlog); | |
| goto end; | |
| } | |
| printf("ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); | |
| fprintf(ficparo,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); | fprintf(ficparo,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); |
| fprintf(ficres,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); | fprintf(ficres,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, 0, weightopt); |
| fprintf(ficlog,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); | fprintf(ficlog,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt); |
| } | } |
| /* ftolpl=6*ftol*1.e5; /\* 6.e-3 make convergences in less than 80 loops for the prevalence limit *\/ */ | /* ftolpl=6*ftol*1.e5; /\* 6.e-3 make convergences in less than 80 loops for the prevalence limit *\/ */ |
| Line 11075 int main(int argc, char *argv[]) | Line 11165 int main(int argc, char *argv[]) |
| ungetc(c,ficpar); | ungetc(c,ficpar); |
| covar=matrix(0,NCOVMAX,1,n); /**< used in readdata */ | covar=matrix(0,NCOVMAX,firstobs,lastobs); /**< used in readdata */ |
| if(nqv>=1)coqvar=matrix(1,nqv,1,n); /**< Fixed quantitative covariate */ | if(nqv>=1)coqvar=matrix(1,nqv,firstobs,lastobs); /**< Fixed quantitative covariate */ |
| if(nqtv>=1)cotqvar=ma3x(1,maxwav,1,nqtv,1,n); /**< Time varying quantitative covariate */ | if(nqtv>=1)cotqvar=ma3x(1,maxwav,1,nqtv,firstobs,lastobs); /**< Time varying quantitative covariate */ |
| if(ntv+nqtv>=1)cotvar=ma3x(1,maxwav,1,ntv+nqtv,1,n); /**< Time varying covariate (dummy and quantitative)*/ | if(ntv+nqtv>=1)cotvar=ma3x(1,maxwav,1,ntv+nqtv,firstobs,lastobs); /**< Time varying covariate (dummy and quantitative)*/ |
| cptcovn=0; /*Number of covariates, i.e. number of '+' in model statement plus one, indepently of n in Vn*/ | cptcovn=0; /*Number of covariates, i.e. number of '+' in model statement plus one, indepently of n in Vn*/ |
| /* v1+v2+v3+v2*v4+v5*age makes cptcovn = 5 | /* v1+v2+v3+v2*v4+v5*age makes cptcovn = 5 |
| v1+v2*age+v2*v3 makes cptcovn = 3 | v1+v2*age+v2*v3 makes cptcovn = 3 |
| Line 11141 int main(int argc, char *argv[]) | Line 11231 int main(int argc, char *argv[]) |
| for(jj=1; jj <=nlstate+ndeath; jj++){ | for(jj=1; jj <=nlstate+ndeath; jj++){ |
| if(jj==i) continue; | if(jj==i) continue; |
| j++; | j++; |
| while((c=getc(ficpar))=='#' && c!= EOF){ | |
| ungetc(c,ficpar); | |
| fgets(line, MAXLINE, ficpar); | |
| numlinepar++; | |
| fputs(line,stdout); | |
| fputs(line,ficparo); | |
| fputs(line,ficlog); | |
| } | |
| ungetc(c,ficpar); | |
| fscanf(ficpar,"%1d%1d",&i1,&j1); | fscanf(ficpar,"%1d%1d",&i1,&j1); |
| if ((i1 != i) || (j1 != jj)){ | if ((i1 != i) || (j1 != jj)){ |
| printf("Error in line parameters number %d, %1d%1d instead of %1d%1d \n \ | printf("Error in line parameters number %d, %1d%1d instead of %1d%1d \n \ |
| Line 11281 Please run with mle=-1 to get a correct | Line 11380 Please run with mle=-1 to get a correct |
| /* Main data | /* Main data |
| */ | */ |
| n= lastobs; | nobs=lastobs-firstobs+1; /* was = lastobs;*/ |
| num=lvector(1,n); | /* num=lvector(1,n); */ |
| moisnais=vector(1,n); | /* moisnais=vector(1,n); */ |
| annais=vector(1,n); | /* annais=vector(1,n); */ |
| moisdc=vector(1,n); | /* moisdc=vector(1,n); */ |
| andc=vector(1,n); | /* andc=vector(1,n); */ |
| weight=vector(1,n); | /* weight=vector(1,n); */ |
| agedc=vector(1,n); | /* agedc=vector(1,n); */ |
| cod=ivector(1,n); | /* cod=ivector(1,n); */ |
| for(i=1;i<=n;i++){ | /* for(i=1;i<=n;i++){ */ |
| num=lvector(firstobs,lastobs); | |
| moisnais=vector(firstobs,lastobs); | |
| annais=vector(firstobs,lastobs); | |
| moisdc=vector(firstobs,lastobs); | |
| andc=vector(firstobs,lastobs); | |
| weight=vector(firstobs,lastobs); | |
| agedc=vector(firstobs,lastobs); | |
| cod=ivector(firstobs,lastobs); | |
| for(i=firstobs;i<=lastobs;i++){ | |
| num[i]=0; | num[i]=0; |
| moisnais[i]=0; | moisnais[i]=0; |
| annais[i]=0; | annais[i]=0; |
| Line 11300 Please run with mle=-1 to get a correct | Line 11408 Please run with mle=-1 to get a correct |
| cod[i]=0; | cod[i]=0; |
| weight[i]=1.0; /* Equal weights, 1 by default */ | weight[i]=1.0; /* Equal weights, 1 by default */ |
| } | } |
| mint=matrix(1,maxwav,1,n); | mint=matrix(1,maxwav,firstobs,lastobs); |
| anint=matrix(1,maxwav,1,n); | anint=matrix(1,maxwav,firstobs,lastobs); |
| s=imatrix(1,maxwav+1,1,n); /* s[i][j] health state for wave i and individual j */ | s=imatrix(1,maxwav+1,firstobs,lastobs); /* s[i][j] health state for wave i and individual j */ |
| tab=ivector(1,NCOVMAX); | tab=ivector(1,NCOVMAX); |
| ncodemax=ivector(1,NCOVMAX); /* Number of code per covariate; if O and 1 only, 2**ncov; V1+V2+V3+V4=>16 */ | ncodemax=ivector(1,NCOVMAX); /* Number of code per covariate; if O and 1 only, 2**ncov; V1+V2+V3+V4=>16 */ |
| ncodemaxwundef=ivector(1,NCOVMAX); /* Number of code per covariate; if - 1 O and 1 only, 2**ncov; V1+V2+V3+V4=>16 */ | ncodemaxwundef=ivector(1,NCOVMAX); /* Number of code per covariate; if - 1 O and 1 only, 2**ncov; V1+V2+V3+V4=>16 */ |
| Line 11404 Please run with mle=-1 to get a correct | Line 11512 Please run with mle=-1 to get a correct |
| agegomp=(int)agemin; | agegomp=(int)agemin; |
| free_vector(moisnais,1,n); | free_vector(moisnais,firstobs,lastobs); |
| free_vector(annais,1,n); | free_vector(annais,firstobs,lastobs); |
| /* free_matrix(mint,1,maxwav,1,n); | /* free_matrix(mint,1,maxwav,1,n); |
| free_matrix(anint,1,maxwav,1,n);*/ | free_matrix(anint,1,maxwav,1,n);*/ |
| /* free_vector(moisdc,1,n); */ | /* free_vector(moisdc,1,n); */ |
| Line 11431 Please run with mle=-1 to get a correct | Line 11539 Please run with mle=-1 to get a correct |
| concatwav(wav, dh, bh, mw, s, agedc, agev, firstpass, lastpass, imx, nlstate, stepm); | concatwav(wav, dh, bh, mw, s, agedc, agev, firstpass, lastpass, imx, nlstate, stepm); |
| /* Concatenates waves */ | /* Concatenates waves */ |
| free_vector(moisdc,1,n); | free_vector(moisdc,firstobs,lastobs); |
| free_vector(andc,1,n); | free_vector(andc,firstobs,lastobs); |
| /* Routine tricode is to calculate cptcoveff (real number of unique covariates) and to associate covariable number and modality */ | /* Routine tricode is to calculate cptcoveff (real number of unique covariates) and to associate covariable number and modality */ |
| nbcode=imatrix(0,NCOVMAX,0,NCOVMAX); | nbcode=imatrix(0,NCOVMAX,0,NCOVMAX); |
| Line 11614 Title=%s <br>Datafile=%s Firstpass=%d La | Line 11722 Title=%s <br>Datafile=%s Firstpass=%d La |
| firstpass, lastpass, stepm, weightopt, model); | firstpass, lastpass, stepm, weightopt, model); |
| fprintf(fichtm,"\n"); | fprintf(fichtm,"\n"); |
| fprintf(fichtm,"<h4>Parameter line 2</h4><ul><li>Tolerance for the convergence of the likelihood: ftol=%f \n<li>Interval for the elementary matrix (in month): stepm=%d",\ | fprintf(fichtm,"<h4>Parameter line 2</h4><ul><li>Tolerance for the convergence of the likelihood: ftol=%g \n<li>Interval for the elementary matrix (in month): stepm=%d",\ |
| ftol, stepm); | ftol, stepm); |
| fprintf(fichtm,"\n<li>Number of fixed dummy covariates: ncovcol=%d ", ncovcol); | fprintf(fichtm,"\n<li>Number of fixed dummy covariates: ncovcol=%d ", ncovcol); |
| ncurrv=1; | ncurrv=1; |
| Line 11622 Title=%s <br>Datafile=%s Firstpass=%d La | Line 11730 Title=%s <br>Datafile=%s Firstpass=%d La |
| fprintf(fichtm,"\n<li> Number of fixed quantitative variables: nqv=%d ", nqv); | fprintf(fichtm,"\n<li> Number of fixed quantitative variables: nqv=%d ", nqv); |
| ncurrv=i; | ncurrv=i; |
| for(i=ncurrv; i <=ncurrv-1+nqv; i++) fprintf(fichtm,"V%d ", i); | for(i=ncurrv; i <=ncurrv-1+nqv; i++) fprintf(fichtm,"V%d ", i); |
| fprintf(fichtm,"\n<li> Number of time varying (wave varying) covariates: ntv=%d ", ntv); | fprintf(fichtm,"\n<li> Number of time varying (wave varying) dummy covariates: ntv=%d ", ntv); |
| ncurrv=i; | ncurrv=i; |
| for(i=ncurrv; i <=ncurrv-1+ntv; i++) fprintf(fichtm,"V%d ", i); | for(i=ncurrv; i <=ncurrv-1+ntv; i++) fprintf(fichtm,"V%d ", i); |
| fprintf(fichtm,"\n<li>Number of quantitative time varying covariates: nqtv=%d ", nqtv); | fprintf(fichtm,"\n<li>Number of time varying quantitative covariates: nqtv=%d ", nqtv); |
| ncurrv=i; | ncurrv=i; |
| for(i=ncurrv; i <=ncurrv-1+nqtv; i++) fprintf(fichtm,"V%d ", i); | for(i=ncurrv; i <=ncurrv-1+nqtv; i++) fprintf(fichtm,"V%d ", i); |
| fprintf(fichtm,"\n<li>Weights column \n<br>Number of alive states: nlstate=%d <br>Number of death states (not really implemented): ndeath=%d \n<li>Number of waves: maxwav=%d \n<li>Parameter for maximization (1), using parameter values (0), for design of parameters and variance-covariance matrix: mle=%d \n<li>Does the weight column be taken into account (1), or not (0): weight=%d</ul>\n", \ | fprintf(fichtm,"\n<li>Weights column \n<br>Number of alive states: nlstate=%d <br>Number of death states (not really implemented): ndeath=%d \n<li>Number of waves: maxwav=%d \n<li>Parameter for maximization (1), using parameter values (0), for design of parameters and variance-covariance matrix: mle=%d \n<li>Does the weight column be taken into account (1), or not (0): weight=%d</ul>\n", \ |
| Line 11657 Interval (in months) between two waves: | Line 11765 Interval (in months) between two waves: |
| for(j=1;j<=NDIM;j++) | for(j=1;j<=NDIM;j++) |
| ximort[i][j]=0.; | ximort[i][j]=0.; |
| /* ximort=gsl_matrix_alloc(1,NDIM,1,NDIM); */ | /* ximort=gsl_matrix_alloc(1,NDIM,1,NDIM); */ |
| cens=ivector(1,n); | cens=ivector(firstobs,lastobs); |
| ageexmed=vector(1,n); | ageexmed=vector(firstobs,lastobs); |
| agecens=vector(1,n); | agecens=vector(firstobs,lastobs); |
| dcwave=ivector(1,n); | dcwave=ivector(firstobs,lastobs); |
| for (i=1; i<=imx; i++){ | for (i=1; i<=imx; i++){ |
| dcwave[i]=-1; | dcwave[i]=-1; |
| Line 11874 Please run with mle=-1 to get a correct | Line 11982 Please run with mle=-1 to get a correct |
| free_vector(lpop,1,AGESUP); | free_vector(lpop,1,AGESUP); |
| free_vector(tpop,1,AGESUP); | free_vector(tpop,1,AGESUP); |
| free_matrix(ximort,1,NDIM,1,NDIM); | free_matrix(ximort,1,NDIM,1,NDIM); |
| free_ivector(cens,1,n); | free_ivector(dcwave,firstobs,lastobs); |
| free_vector(agecens,1,n); | free_vector(agecens,firstobs,lastobs); |
| free_ivector(dcwave,1,n); | free_vector(ageexmed,firstobs,lastobs); |
| free_ivector(cens,firstobs,lastobs); | |
| #ifdef GSL | #ifdef GSL |
| #endif | #endif |
| } /* Endof if mle==-3 mortality only */ | } /* Endof if mle==-3 mortality only */ |
| Line 12282 Please run with mle=-1 to get a correct | Line 12391 Please run with mle=-1 to get a correct |
| /* free_imatrix(dh,1,lastpass-firstpass+2,1,imx); */ | /* free_imatrix(dh,1,lastpass-firstpass+2,1,imx); */ |
| /* free_imatrix(bh,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_imatrix(mw,1,lastpass-firstpass+2,1,imx); */ |
| free_lvector(num,1,n); | free_lvector(num,firstobs,lastobs); |
| free_vector(agedc,1,n); | free_vector(agedc,firstobs,lastobs); |
| /*free_matrix(covar,0,NCOVMAX,1,n);*/ | /*free_matrix(covar,0,NCOVMAX,1,n);*/ |
| /*free_matrix(covar,1,NCOVMAX,1,n);*/ | /*free_matrix(covar,1,NCOVMAX,1,n);*/ |
| fclose(ficparo); | fclose(ficparo); |
| Line 12540 Please run with mle=-1 to get a correct | Line 12649 Please run with mle=-1 to get a correct |
| if(vpopbased==1) | if(vpopbased==1) |
| fprintf(ficrest,"the age specific prevalence observed (cross-sectionally) in the population i.e cross-sectionally\n in each health state (popbased=1) (mobilav=%d)\n",mobilav); | fprintf(ficrest,"the age specific prevalence observed (cross-sectionally) in the population i.e cross-sectionally\n in each health state (popbased=1) (mobilav=%d)\n",mobilav); |
| else | else |
| fprintf(ficrest,"the age specific period (stable) prevalences in each health state \n"); | fprintf(ficrest,"the age specific forward period (stable) prevalences in each health state \n"); |
| fprintf(ficrest,"# Age popbased mobilav e.. (std) "); | fprintf(ficrest,"# Age popbased mobilav e.. (std) "); |
| for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i); | for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i); |
| fprintf(ficrest,"\n"); | fprintf(ficrest,"\n"); |
| /* printf("Which p?\n"); for(i=1;i<=npar;i++)printf("p[i=%d]=%lf,",i,p[i]);printf("\n"); */ | /* printf("Which p?\n"); for(i=1;i<=npar;i++)printf("p[i=%d]=%lf,",i,p[i]);printf("\n"); */ |
| printf("Computing age specific period (stable) prevalences in each health state \n"); | printf("Computing age specific forward period (stable) prevalences in each health state \n"); |
| fprintf(ficlog,"Computing age specific period (stable) prevalences in each health state \n"); | fprintf(ficlog,"Computing age specific forward period (stable) prevalences in each health state \n"); |
| for(age=bage; age <=fage ;age++){ | for(age=bage; age <=fage ;age++){ |
| prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, &ncvyear, k, nres); /*ZZ Is it the correct prevalim */ | prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, &ncvyear, k, nres); /*ZZ Is it the correct prevalim */ |
| if (vpopbased==1) { | if (vpopbased==1) { |
| Line 12593 Please run with mle=-1 to get a correct | Line 12702 Please run with mle=-1 to get a correct |
| printf("done State-specific expectancies\n");fflush(stdout); | printf("done State-specific expectancies\n");fflush(stdout); |
| fprintf(ficlog,"done State-specific expectancies\n");fflush(ficlog); | fprintf(ficlog,"done State-specific expectancies\n");fflush(ficlog); |
| /* variance-covariance of period prevalence*/ | /* variance-covariance of forward period prevalence*/ |
| varprlim(fileresu, nresult, mobaverage, mobilavproj, bage, fage, prlim, &ncvyear, ftolpl, p, matcov, delti, stepm, cptcoveff); | varprlim(fileresu, nresult, mobaverage, mobilavproj, bage, fage, prlim, &ncvyear, ftolpl, p, matcov, delti, stepm, cptcoveff); |
| free_vector(weight,1,n); | free_vector(weight,firstobs,lastobs); |
| free_imatrix(Tvard,1,NCOVMAX,1,2); | free_imatrix(Tvard,1,NCOVMAX,1,2); |
| free_imatrix(s,1,maxwav+1,1,n); | free_imatrix(s,1,maxwav+1,firstobs,lastobs); |
| free_matrix(anint,1,maxwav,1,n); | free_matrix(anint,1,maxwav,firstobs,lastobs); |
| free_matrix(mint,1,maxwav,1,n); | free_matrix(mint,1,maxwav,firstobs,lastobs); |
| free_ivector(cod,1,n); | free_ivector(cod,firstobs,lastobs); |
| free_ivector(tab,1,NCOVMAX); | free_ivector(tab,1,NCOVMAX); |
| fclose(ficresstdeij); | fclose(ficresstdeij); |
| fclose(ficrescveij); | fclose(ficrescveij); |
| Line 12622 Please run with mle=-1 to get a correct | Line 12731 Please run with mle=-1 to get a correct |
| free_matrix(oldms, 1,nlstate+ndeath,1,nlstate+ndeath); | free_matrix(oldms, 1,nlstate+ndeath,1,nlstate+ndeath); |
| free_matrix(newms, 1,nlstate+ndeath,1,nlstate+ndeath); | free_matrix(newms, 1,nlstate+ndeath,1,nlstate+ndeath); |
| free_matrix(savms, 1,nlstate+ndeath,1,nlstate+ndeath); | free_matrix(savms, 1,nlstate+ndeath,1,nlstate+ndeath); |
| if(ntv+nqtv>=1)free_ma3x(cotvar,1,maxwav,1,ntv+nqtv,1,n); | if(ntv+nqtv>=1)free_ma3x(cotvar,1,maxwav,1,ntv+nqtv,firstobs,lastobs); |
| if(nqtv>=1)free_ma3x(cotqvar,1,maxwav,1,nqtv,1,n); | if(nqtv>=1)free_ma3x(cotqvar,1,maxwav,1,nqtv,firstobs,lastobs); |
| if(nqv>=1)free_matrix(coqvar,1,nqv,1,n); | if(nqv>=1)free_matrix(coqvar,1,nqv,firstobs,lastobs); |
| free_matrix(covar,0,NCOVMAX,1,n); | free_matrix(covar,0,NCOVMAX,firstobs,lastobs); |
| free_matrix(matcov,1,npar,1,npar); | free_matrix(matcov,1,npar,1,npar); |
| free_matrix(hess,1,npar,1,npar); | free_matrix(hess,1,npar,1,npar); |
| /*free_vector(delti,1,npar);*/ | /*free_vector(delti,1,npar);*/ |
| Line 12745 Please run with mle=-1 to get a correct | Line 12854 Please run with mle=-1 to get a correct |
| sprintf(plotcmd,"%s %s",pplotcmd, optionfilegnuplot); | sprintf(plotcmd,"%s %s",pplotcmd, optionfilegnuplot); |
| printf("Starting graphs with: '%s'\n",plotcmd);fflush(stdout); | printf("Starting graphs with: '%s'\n",plotcmd);fflush(stdout); |
| strcpy(pplotcmd,plotcmd); | |
| if((outcmd=system(plotcmd)) != 0){ | if((outcmd=system(plotcmd)) != 0){ |
| printf("gnuplot command might not be in your path: '%s', err=%d\n", plotcmd, outcmd); | printf("Error in gnuplot, command might not be in your path: '%s', err=%d\n", plotcmd, outcmd); |
| printf("\n Trying if gnuplot resides on the same directory that IMaCh\n"); | printf("\n Trying if gnuplot resides on the same directory that IMaCh\n"); |
| sprintf(plotcmd,"%sgnuplot %s", pathimach, optionfilegnuplot); | sprintf(plotcmd,"%sgnuplot %s", pathimach, optionfilegnuplot); |
| if((outcmd=system(plotcmd)) != 0) | if((outcmd=system(plotcmd)) != 0){ |
| printf("\n Still a problem with gnuplot command %s, err=%d\n", plotcmd, outcmd); | printf("\n Still a problem with gnuplot command %s, err=%d\n", plotcmd, outcmd); |
| strcpy(plotcmd,pplotcmd); | |
| } | |
| } | } |
| printf(" Successful, please wait..."); | printf(" Successful, please wait..."); |
| while (z[0] != 'q') { | while (z[0] != 'q') { |