Imach, Version 0.71e
Title=%s
Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=%s
Total number of observations=%d
@@ -1910,8 +1919,9 @@ Interval (in months) between two waves:
- Variances of life expectancies by age and initial health status: v%s
- Health expectancies with their variances: t%s
- Standard deviation of stationary prevalences: vpl%s
- - Prevalences and population forecasting: f%s
-
",title,datafile,firstpass,lastpass,stepm, weightopt,model,imx,jmin,jmax,jmean,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres);
+ - Prevalences forecasting: f%s
+ - Population forecasting (if popforecast=1): pop%s
+
",title,datafile,firstpass,lastpass,stepm, weightopt,model,imx,jmin,jmax,jmean,fileres,fileres,rfileres,rfileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres,fileres);
fprintf(fichtm," - Graphs
");
@@ -1953,7 +1963,7 @@ fclose(fichtm);
}
/******************* Gnuplot file **************/
-void printinggnuplot(char fileres[],char optionfilefiname[],char optionfile[],char optionfilegnuplot[], double agemin, double agemax, double fage , char pathc[], double p[]){
+void printinggnuplot(char fileres[],char optionfilefiname[],char optionfile[],char optionfilegnuplot[], double agemin, double agemaxpar, double fage , char pathc[], double p[]){
int m,cpt,k1,i,k,j,jk,k2,k3,ij,l;
@@ -2048,7 +2058,7 @@ fprintf(ficgp,"\nset out \"v%s%d%d.gif\"
for (k1=1; k1<= m ; k1 ++) {
for (cpt=1; cpt=(agemin-((int)calagedate %12)/12.); agedeb--){
+ nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm);
+ nhstepm = nhstepm/hstepm;
+
+ p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ oldm=oldms;savm=savms;
+ hpxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k);
+
+ for (h=0; h<=nhstepm; h++){
+ if (h==(int) (calagedate+YEARM*cpt)) {
+ fprintf(ficresf,"\n %.f ",agedeb+h*hstepm/YEARM*stepm);
+ }
+ for(j=1; j<=nlstate+ndeath;j++) {
+ kk1=0.;kk2=0;
+ for(i=1; i<=nlstate;i++) {
+ if (mobilav==1)
+ kk1=kk1+p3mat[i][j][h]*mobaverage[(int)agedeb+1][i][cptcod];
+ else {
+ kk1=kk1+p3mat[i][j][h]*probs[(int)(agedeb+1)][i][cptcod];
+ }
+
+ }
+ if (h==(int)(calagedate+12*cpt)){
+ fprintf(ficresf," %.3f", kk1);
+
+ }
+ }
+ }
+ free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ }
+ }
+ }
+ }
+
+ if (mobilav==1) free_ma3x(mobaverage,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+
+ fclose(ficresf);
+}
+/************** Forecasting ******************/
+populforecast(char fileres[], double anpyram,double mpyram,double jpyram,double agemin, double agemax,double dateprev1, double dateprev2, int mobilav, double agedeb, double fage, int popforecast, char popfile[], double anpyram1,double p[], int i2){
+
+ int cpt, stepsize, hstepm, nhstepm, j,k,c, cptcod, i,h;
+ int *popage;
+ double calagedate, agelim, kk1, kk2, yp,yp1,yp2,jprojmean,mprojmean,anprojmean;
+ double *popeffectif,*popcount;
+ double ***p3mat,***tabpop,***tabpopprev;
+ char filerespop[FILENAMELENGTH];
+
+ tabpop= ma3x(1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+ tabpopprev= ma3x(1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+ agelim=AGESUP;
+ calagedate=(anpyram+mpyram/12.+jpyram/365.-dateintmean)*YEARM;
+
+ prevalence(agemin, agemax, s, agev, nlstate, imx,Tvar,nbcode, ncodemax,mint,anint,dateprev1,dateprev2, calagedate);
+
+
+ strcpy(filerespop,"pop");
+ strcat(filerespop,fileres);
+ if((ficrespop=fopen(filerespop,"w"))==NULL) {
+ printf("Problem with forecast resultfile: %s\n", filerespop);
+ }
+ printf("Computing forecasting: result on file '%s' \n", filerespop);
+
+ if (cptcoveff==0) ncodemax[cptcoveff]=1;
+
+ if (mobilav==1) {
+ mobaverage= ma3x(1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+ movingaverage(agedeb, fage, agemin, mobaverage);
+ }
+
+ stepsize=(int) (stepm+YEARM-1)/YEARM;
+ if (stepm<=12) stepsize=1;
+
+ agelim=AGESUP;
+
+ hstepm=1;
+ hstepm=hstepm/stepm;
+
+ if (popforecast==1) {
+ if((ficpop=fopen(popfile,"r"))==NULL) {
+ printf("Problem with population file : %s\n",popfile);exit(0);
+ }
+ popage=ivector(0,AGESUP);
+ popeffectif=vector(0,AGESUP);
+ popcount=vector(0,AGESUP);
+
+ i=1;
+ while ((c=fscanf(ficpop,"%d %lf\n",&popage[i],&popcount[i])) != EOF) i=i+1;
+
+ imx=i;
+ for (i=1; i=(agemin-((int)calagedate %12)/12.); agedeb--){
+ nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm);
+ nhstepm = nhstepm/hstepm;
+
+ p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ oldm=oldms;savm=savms;
+ hpxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k);
+
+ for (h=0; h<=nhstepm; h++){
+ if (h==(int) (calagedate+YEARM*cpt)) {
+ fprintf(ficrespop,"\n %3.f ",agedeb+h*hstepm/YEARM*stepm);
+ }
+ for(j=1; j<=nlstate+ndeath;j++) {
+ kk1=0.;kk2=0;
+ for(i=1; i<=nlstate;i++) {
+ if (mobilav==1)
+ kk1=kk1+p3mat[i][j][h]*mobaverage[(int)agedeb+1][i][cptcod];
+ else {
+ kk1=kk1+p3mat[i][j][h]*probs[(int)(agedeb+1)][i][cptcod];
+ }
+ }
+ if (h==(int)(calagedate+12*cpt)){
+ tabpop[(int)(agedeb)][j][cptcod]=kk1;
+ /*fprintf(ficrespop," %.3f", kk1);
+ if (popforecast==1) fprintf(ficrespop," [%.f]", kk1*popeffectif[(int)agedeb+1]);*/
+ }
+ }
+ for(i=1; i<=nlstate;i++){
+ kk1=0.;
+ for(j=1; j<=nlstate;j++){
+ kk1= kk1+tabpop[(int)(agedeb)][j][cptcod];
+ }
+ tabpopprev[(int)(agedeb)][i][cptcod]=tabpop[(int)(agedeb)][i][cptcod]/kk1*popeffectif[(int)(agedeb+(calagedate+12*cpt)*hstepm/YEARM*stepm-1)];
+ }
+
+ if (h==(int)(calagedate+12*cpt)) for(j=1; j<=nlstate;j++)
+ fprintf(ficrespop," %15.2f",tabpopprev[(int)(agedeb+1)][j][cptcod]);
+ }
+ free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ }
+ }
+
+ /******/
+
+ for (cpt=1; cpt<=(anpyram1-anpyram);cpt++) {
+ fprintf(ficrespop,"\n\n# Forecasting at date %.lf/%.lf/%.lf ",jpyram,mpyram,anpyram+cpt);
+ for (agedeb=(fage-((int)calagedate %12/12.)); agedeb>=(agemin-((int)calagedate %12)/12.); agedeb--){
+ nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm);
+ nhstepm = nhstepm/hstepm;
+
+ p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ oldm=oldms;savm=savms;
+ hpxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k);
+ for (h=0; h<=nhstepm; h++){
+ if (h==(int) (calagedate+YEARM*cpt)) {
+ fprintf(ficresf,"\n %3.f ",agedeb+h*hstepm/YEARM*stepm);
+ }
+ for(j=1; j<=nlstate+ndeath;j++) {
+ kk1=0.;kk2=0;
+ for(i=1; i<=nlstate;i++) {
+ kk1=kk1+p3mat[i][j][h]*tabpopprev[(int)agedeb+1][i][cptcod];
+ }
+ if (h==(int)(calagedate+12*cpt)) fprintf(ficresf," %15.2f", kk1);
+ }
+ }
+ free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
+ }
+ }
+ }
+ }
+
+ if (mobilav==1) free_ma3x(mobaverage,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+
+ if (popforecast==1) {
+ free_ivector(popage,0,AGESUP);
+ free_vector(popeffectif,0,AGESUP);
+ free_vector(popcount,0,AGESUP);
+ }
+ free_ma3x(tabpop,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+ free_ma3x(tabpopprev,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
+ fclose(ficrespop);
+}
+
/***********************************************/
/**************** Main Program *****************/
/***********************************************/
@@ -2155,7 +2421,7 @@ int main(int argc, char *argv[])
int i,j, k, n=MAXN,iter,m,size,cptcode, cptcod;
double agedeb, agefin,hf;
- double agemin=1.e20, agemax=-1.e20;
+ double agemin=1.e20, agemaxpar=-1.e20, agemax=-1.e20;
double fret;
double **xi,tmp,delta;
@@ -2168,7 +2434,7 @@ int main(int argc, char *argv[])
char optionfile[FILENAMELENGTH], datafile[FILENAMELENGTH], filerespl[FILENAMELENGTH];
char optionfilext[10], optionfilefiname[FILENAMELENGTH], optionfilehtm[FILENAMELENGTH], optionfilegnuplot[FILENAMELENGTH], plotcmd[FILENAMELENGTH];
- char fileres[FILENAMELENGTH], filerespij[FILENAMELENGTH], filereso[FILENAMELENGTH], fileresf[FILENAMELENGTH];
+ char fileres[FILENAMELENGTH], filerespij[FILENAMELENGTH], filereso[FILENAMELENGTH], rfileres[FILENAMELENGTH];
char filerest[FILENAMELENGTH];
char fileregp[FILENAMELENGTH];
@@ -2182,8 +2448,7 @@ int main(int argc, char *argv[])
int jnais,jdc,jint4,jint1,jint2,jint3,**outcome,**adl,*tab;
int mobilav=0,popforecast=0;
int hstepm, nhstepm;
- int *popage;/*boolprev=0 if date and zero if wave*/
- double jprev1, mprev1,anprev1,jprev2, mprev2,anprev2;
+ double jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,jpyram, mpyram,anpyram,jpyram1, mpyram1,anpyram1;
double bage, fage, age, agelim, agebase;
double ftolpl=FTOL;
@@ -2198,11 +2463,10 @@ int main(int argc, char *argv[])
double **varpl; /* Variances of prevalence limits by age */
double *epj, vepp;
double kk1, kk2;
- double *popeffectif,*popcount;
- double dateprev1, dateprev2,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,jprojmean,mprojmean,anprojmean, calagedate;
- double yp,yp1,yp2;
+ double dateprev1, dateprev2,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2;
+
- char version[80]="Imach version 0.7, February 2002, INED-EUROREVES ";
+ char version[80]="Imach version 0.71e, March 2002, INED-EUROREVES ";
char *alph[]={"a","a","b","c","d","e"}, str[4];
@@ -2370,12 +2634,17 @@ while((c=getc(ficpar))=='#' && c!= EOF){
printf("\n");
- /*-------- data file ----------*/
- if((ficres =fopen(fileres,"w"))==NULL) {
- printf("Problem with resultfile: %s\n", fileres);goto end;
+ /*-------- Rewriting paramater file ----------*/
+ strcpy(rfileres,"r"); /* "Rparameterfile */
+ strcat(rfileres,optionfilefiname); /* Parameter file first name*/
+ strcat(rfileres,"."); /* */
+ strcat(rfileres,optionfilext); /* Other files have txt extension */
+ if((ficres =fopen(rfileres,"w"))==NULL) {
+ printf("Problem writing new parameter file: %s\n", fileres);goto end;
}
fprintf(ficres,"#%s\n",version);
+ /*-------- data file ----------*/
if((fic=fopen(datafile,"r"))==NULL) {
printf("Problem with datafile: %s\n", datafile);goto end;
}
@@ -2648,7 +2917,7 @@ printf("Total number of individuals= %d,
}
/*--------- results files --------------*/
- fprintf(ficres,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncov=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\nmodel=%s\n", title, datafile, lastobs, firstpass,lastpass,ftol, stepm, ncov, nlstate, ndeath, maxwav, mle,weightopt,model);
+ fprintf(ficres,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncov=%d nlstate=%d ndeath=%d maxwav=%d mle= 0 weight=%d\nmodel=%s\n", title, datafile, lastobs, firstpass,lastpass,ftol, stepm, ncov, nlstate, ndeath, maxwav, weightopt,model);
jk=1;
@@ -2720,16 +2989,16 @@ printf("Total number of individuals= %d,
}
ungetc(c,ficpar);
- fscanf(ficpar,"agemin=%lf agemax=%lf bage=%lf fage=%lf\n",&agemin,&agemax, &bage, &fage);
+ fscanf(ficpar,"agemin=%lf agemax=%lf bage=%lf fage=%lf\n",&agemin,&agemaxpar, &bage, &fage);
if (fage <= 2) {
bage = agemin;
- fage = agemax;
+ fage = agemaxpar;
}
fprintf(ficres,"# agemin agemax for life expectancy, bage fage (if mle==0 ie no data nor Max likelihood).\n");
- fprintf(ficres,"agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f\n",agemin,agemax,bage,fage);
- fprintf(ficparo,"agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f\n",agemin,agemax,bage,fage);
+ fprintf(ficres,"agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f\n",agemin,agemaxpar,bage,fage);
+ fprintf(ficparo,"agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f\n",agemin,agemaxpar,bage,fage);
while((c=getc(ficpar))=='#' && c!= EOF){
ungetc(c,ficpar);
@@ -2756,9 +3025,9 @@ printf("Total number of individuals= %d,
dateprev2=anprev2+mprev2/12.+jprev2/365.;
fscanf(ficpar,"pop_based=%d\n",&popbased);
- fprintf(ficparo,"pop_based=%d\n",popbased);
- fprintf(ficres,"pop_based=%d\n",popbased);
-
+ fprintf(ficparo,"pop_based=%d\n",popbased);
+ fprintf(ficres,"pop_based=%d\n",popbased);
+
while((c=getc(ficpar))=='#' && c!= EOF){
ungetc(c,ficpar);
fgets(line, MAXLINE, ficpar);
@@ -2766,14 +3035,28 @@ printf("Total number of individuals= %d,
fputs(line,ficparo);
}
ungetc(c,ficpar);
- fscanf(ficpar,"popforecast=%d popfile=%s starting-proj-date=%lf/%lf/%lf final-proj-date=%lf/%lf/%lf mov_average=%d\n",&popforecast,popfile,&jproj1,&mproj1,&anproj1,&jproj2,&mproj2,&anproj2,&mobilav);
-fprintf(ficparo,"popforecast=%d popfile=%s starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mov_average=%d\n",popforecast,popfile,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilav);
-fprintf(ficres,"popforecast=%d popfile=%s starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mov_average=%d\n",popforecast,popfile,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilav);
+
+ fscanf(ficpar,"starting-proj-date=%lf/%lf/%lf final-proj-date=%lf/%lf/%lf mov_average=%d\n",&jproj1,&mproj1,&anproj1,&jproj2,&mproj2,&anproj2,&mobilav);
+fprintf(ficparo,"starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mov_average=%d\n",jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilav);
+fprintf(ficres,"starting-proj-date=%.lf/%.lf/%.lf final-proj-date=%.lf/%.lf/%.lf mov_average=%d\n",jproj1,mproj1,anproj1,jproj2,mproj2,anproj2,mobilav);
+
+
+while((c=getc(ficpar))=='#' && c!= EOF){
+ ungetc(c,ficpar);
+ fgets(line, MAXLINE, ficpar);
+ puts(line);
+ fputs(line,ficparo);
+ }
+ ungetc(c,ficpar);
+
+ fscanf(ficpar,"popforecast=%d popfile=%s popfiledate=%lf/%lf/%lf last-popfiledate=%lf/%lf/%lf\n",&popforecast,popfile,&jpyram,&mpyram,&anpyram,&jpyram1,&mpyram1,&anpyram1);
+ fprintf(ficparo,"popforecast=%d popfile=%s popfiledate=%.lf/%.lf/%.lf last-popfiledate=%.lf/%.lf/%.lf\n",popforecast,popfile,jpyram,mpyram,anpyram,jpyram1,mpyram1,anpyram1);
+ fprintf(ficres,"popforecast=%d popfile=%s popfiledate=%.lf/%.lf/%.lf last-popfiledate=%.lf/%.lf/%.lf\n",popforecast,popfile,jpyram,mpyram,anpyram,jpyram1,mpyram1,anpyram1);
freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvar,nbcode, ncodemax,mint,anint,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);
/*------------ gnuplot -------------*/
- printinggnuplot(fileres,optionfilefiname,optionfile,optionfilegnuplot, agemin,agemax,fage, pathc,p);
+ printinggnuplot(fileres,optionfilefiname,optionfile,optionfilegnuplot, agemin,agemaxpar,fage, pathc,p);
/*------------ free_vector -------------*/
chdir(path);
@@ -2786,22 +3069,10 @@ fprintf(ficres,"popforecast=%d popfile=%
/*free_matrix(covar,1,NCOVMAX,1,n);*/
fclose(ficparo);
fclose(ficres);
-
- /* Reads comments: lines beginning with '#' */
- while((c=getc(ficpar))=='#' && c!= EOF){
- ungetc(c,ficpar);
- fgets(line, MAXLINE, ficpar);
- puts(line);
- fputs(line,ficparo);
- }
- ungetc(c,ficpar);
-
- fscanf(ficpar,"agemin=%lf agemax=%lf bage=%lf fage=%lf\n",&agemin,&agemax, &bage, &fage);
- printf("agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f\n",agemin,agemax, bage, fage);
- fprintf(ficparo,"agemin=%.0f agemax=%.0f bage=%.0f fage=%.0f\n",agemin,agemax,bage,fage);
+
/*--------- index.htm --------*/
- printinghtml(fileres,title,datafile, firstpass, lastpass, stepm, weightopt,model,imx,jmin,jmax,jmean,optionfile,optionfilehtm);
+ printinghtml(fileres,title,datafile, firstpass, lastpass, stepm, weightopt,model,imx,jmin,jmax,jmean,optionfile,optionfilehtm,rfileres);
/*--------------- Prevalence limit --------------*/
@@ -2825,7 +3096,7 @@ fprintf(ficres,"popforecast=%d popfile=%
oldm=oldms; newm=newms; savm=savms; /* Keeps fixed addresses to free */
k=0;
agebase=agemin;
- agelim=agemax;
+ agelim=agemaxpar;
ftolpl=1.e-10;
i1=cptcoveff;
if (cptcovn < 1){i1=1;}
@@ -2902,139 +3173,20 @@ fprintf(ficres,"popforecast=%d popfile=%
fclose(ficrespij);
- if(stepm == 1) {
- /*---------- Forecasting ------------------*/
- calagedate=(anproj1+mproj1/12.+jproj1/365.-dateintmean)*YEARM;
-
-
- prevalence(agemin, agemax, s, agev, nlstate, imx,Tvar,nbcode, ncodemax,mint,anint,dateprev1,dateprev2, calagedate);
-
-
- strcpy(fileresf,"f");
- strcat(fileresf,fileres);
- if((ficresf=fopen(fileresf,"w"))==NULL) {
- printf("Problem with forecast resultfile: %s\n", fileresf);goto end;
- }
- printf("Computing forecasting: result on file '%s' \n", fileresf);
-
- free_matrix(mint,1,maxwav,1,n);
- free_matrix(anint,1,maxwav,1,n);
-
- if (cptcoveff==0) ncodemax[cptcoveff]=1;
-
- if (mobilav==1) {
- mobaverage= ma3x(1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
- movingaverage(agedeb, fage, agemin, mobaverage);
- }
-
- stepsize=(int) (stepm+YEARM-1)/YEARM;
- if (stepm<=12) stepsize=1;
-
- agelim=AGESUP;
-
- hstepm=1;
- hstepm=hstepm/stepm;
- yp1=modf(dateintmean,&yp);
- anprojmean=yp;
- yp2=modf((yp1*12),&yp);
- mprojmean=yp;
- yp1=modf((yp2*30.5),&yp);
- jprojmean=yp;
- if(jprojmean==0) jprojmean=1;
- if(mprojmean==0) jprojmean=1;
-
- fprintf(ficresf,"# Estimated date of observed prevalence: %.lf/%.lf/%.lf ",jprojmean,mprojmean,anprojmean);
-
- if (popforecast==1) {
- if((ficpop=fopen(popfile,"r"))==NULL) {
- printf("Problem with population file : %s\n",popfile);goto end;
- }
- popage=ivector(0,AGESUP);
- popeffectif=vector(0,AGESUP);
- popcount=vector(0,AGESUP);
-
- i=1;
- while ((c=fscanf(ficpop,"%d %lf\n",&popage[i],&popcount[i])) != EOF)
- {
- i=i+1;
- }
- imx=i;
-
- for (i=1; i=(agemin-((int)calagedate %12)/12.); agedeb--){
- nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm);
- nhstepm = nhstepm/hstepm;
-
- p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
- oldm=oldms;savm=savms;
- hpxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k);
-
- for (h=0; h<=nhstepm; h++){
- if (h==(int) (calagedate+YEARM*cpt)) {
- fprintf(ficresf,"\n %.f ",agedeb+h*hstepm/YEARM*stepm);
- }
- for(j=1; j<=nlstate+ndeath;j++) {
- kk1=0.;kk2=0;
- for(i=1; i<=nlstate;i++) {
- if (mobilav==1)
- kk1=kk1+p3mat[i][j][h]*mobaverage[(int)agedeb+1][i][cptcod];
- else {
- kk1=kk1+p3mat[i][j][h]*probs[(int)(agedeb+1)][i][cptcod];
- }
- if (popforecast==1) kk2=kk1*popeffectif[(int)agedeb];
- }
- if (h==(int)(calagedate+12*cpt)){
- fprintf(ficresf," %.3f", kk1);
-
- if (popforecast==1) fprintf(ficresf," [%.f]", kk2);
- }
- }
- }
- free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
- }
- }
- }
- }
-
- if (mobilav==1) free_ma3x(mobaverage,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
-
- if (popforecast==1) {
- free_ivector(popage,0,AGESUP);
- free_vector(popeffectif,0,AGESUP);
- free_vector(popcount,0,AGESUP);
- }
-
- free_imatrix(s,1,maxwav+1,1,n);
- free_vector(weight,1,n);
- fclose(ficresf);
- }
-
- /* End forecasting */
+ /*---------- Forecasting ------------------*/
+ if((stepm == 1) && (strcmp(model,".")==0)){
+ prevforecast(fileres, anproj1,mproj1,jproj1, agemin,agemax, dateprev1, dateprev2,mobilav, agedeb, fage, popforecast, popfile, anproj2,p, i1);
+ if (popforecast==1) populforecast(fileres, anpyram,mpyram,jpyram, agemin,agemax, dateprev1, dateprev2,mobilav, agedeb, fage, popforecast, popfile, anpyram1,p, i1);
+ free_matrix(mint,1,maxwav,1,n);
+ free_matrix(anint,1,maxwav,1,n); free_imatrix(s,1,maxwav+1,1,n);
+ free_vector(weight,1,n);}
else{
erreur=108;
- printf("Error %d!! You can only forecast the prevalences if the optimization\n has been performed with stepm = 1 (month) instead of %d\n", erreur, stepm);
+ printf("Warning %d!! You can only forecast the prevalences if the optimization\n has been performed with stepm = 1 (month) instead of %d or model=. instead of '%s'\n", erreur, stepm, model);
}
-
+
/*---------- Health expectancies and variances ------------*/
strcpy(filerest,"t");
@@ -3091,8 +3243,6 @@ fprintf(ficres,"popforecast=%d popfile=%
for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i);
fprintf(ficrest,"\n");
- hf=1;
- if (stepm >= YEARM) hf=stepm/YEARM;
epj=vector(1,nlstate+1);
for(age=bage; age <=fage ;age++){
prevalim(prlim, nlstate, p, age, oldm, savm,ftolpl,k);
@@ -3101,38 +3251,34 @@ fprintf(ficres,"popforecast=%d popfile=%
prlim[i][i]=probs[(int)age][i][k];
}
- fprintf(ficrest," %.0f",age);
+ fprintf(ficrest," %4.0f",age);
for(j=1, epj[nlstate+1]=0.;j <=nlstate;j++){
for(i=1, epj[j]=0.;i <=nlstate;i++) {
- epj[j] += prlim[i][i]*hf*eij[i][j][(int)age];
+ epj[j] += prlim[i][i]*eij[i][j][(int)age];
}
epj[nlstate+1] +=epj[j];
}
for(i=1, vepp=0.;i <=nlstate;i++)
for(j=1;j <=nlstate;j++)
vepp += vareij[i][j][(int)age];
- fprintf(ficrest," %.2f (%.2f)", epj[nlstate+1],hf*sqrt(vepp));
+ fprintf(ficrest," %7.2f (%7.2f)", epj[nlstate+1],sqrt(vepp));
for(j=1;j <=nlstate;j++){
- fprintf(ficrest," %.2f (%.2f)", epj[j],hf*sqrt(vareij[j][j][(int)age]));
+ fprintf(ficrest," %7.2f (%7.2f)", epj[j],sqrt(vareij[j][j][(int)age]));
}
fprintf(ficrest,"\n");
}
}
}
-
-
-
- fclose(ficreseij);
- fclose(ficresvij);
+ fclose(ficreseij);
+ fclose(ficresvij);
fclose(ficrest);
fclose(ficpar);
free_vector(epj,1,nlstate+1);
- /* scanf("%d ",i); */
/*------- Variance limit prevalence------*/
-strcpy(fileresvpl,"vpl");
+ strcpy(fileresvpl,"vpl");
strcat(fileresvpl,fileres);
if((ficresvpl=fopen(fileresvpl,"w"))==NULL) {
printf("Problem with variance prev lim resultfile: %s\n", fileresvpl);
@@ -3140,19 +3286,19 @@ strcpy(fileresvpl,"vpl");
}
printf("Computing Variance-covariance of Prevalence limit: file '%s' \n", fileresvpl);
- k=0;
- for(cptcov=1;cptcov<=i1;cptcov++){
- for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){
- k=k+1;
- fprintf(ficresvpl,"\n#****** ");
- for(j=1;j<=cptcoveff;j++)
- fprintf(ficresvpl,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtab[k][j]]);
- fprintf(ficresvpl,"******\n");
-
- varpl=matrix(1,nlstate,(int) bage, (int) fage);
- oldm=oldms;savm=savms;
+ k=0;
+ for(cptcov=1;cptcov<=i1;cptcov++){
+ for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){
+ k=k+1;
+ fprintf(ficresvpl,"\n#****** ");
+ for(j=1;j<=cptcoveff;j++)
+ fprintf(ficresvpl,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtab[k][j]]);
+ fprintf(ficresvpl,"******\n");
+
+ varpl=matrix(1,nlstate,(int) bage, (int) fage);
+ oldm=oldms;savm=savms;
varprevlim(fileres, varpl, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k);
- }
+ }
}
fclose(ficresvpl);