/* $Id$
$State$
$Log$
+ Revision 1.350 2023/04/24 11:38:06 brouard
+ *** empty log message ***
+
Revision 1.349 2023/01/31 09:19:37 brouard
Summary: Improvements in models with age*Vn*Vm
return s;
}
+char *trimbtab(char *out, char *in)
+{ /* Trim blanks or tabs in line but keeps first blanks if line starts with blanks */
+ char *s;
+ s=out;
+ while (*in != '\0'){
+ while( (*in == ' ' || *in == '\t')){ /* && *(in+1) != '\0'){*/
+ in++;
+ }
+ *out++ = *in++;
+ }
+ *out='\0';
+ return s;
+}
+
/* char *substrchaine(char *out, char *in, char *chain) */
/* { */
/* /\* Substract chain 'chain' from 'in', return and output 'out' *\/ */
fprintf(ficgp," u %d:(",ioffset);
kl=0;
strcpy(gplotcondition,"(");
- for (k=1; k<=cptcoveff; k++){ /* For each covariate writing the chain of conditions */
+ /* for (k=1; k<=cptcoveff; k++){ /\* For each covariate writing the chain of conditions *\/ */
/* lv= decodtabm(k1,k,cptcoveff); /\* Should be the covariate value corresponding to combination k1 and covariate k *\/ */
- lv=codtabm(k1,TnsdVar[Tvaraff[k]]);
+ for (k=1; k<=cptcovs; k++){ /* For each covariate k get corresponding value lv for combination k1 */
+ /* lv=codtabm(k1,TnsdVar[Tvaraff[k]]); */
+ lv=Tvresult[nres][k];
+ vlv=TinvDoQresult[nres][Tvresult[nres][k]];
/* decodtabm(1,1,4) = 1 because h=1 k= (1) 1 1 1 */
/* decodtabm(1,2,4) = 1 because h=1 k= 1 (1) 1 1 */
/* decodtabm(13,3,4)= 2 because h=13 k= 1 1 (2) 2 */
/* vlv= nbcode[Tvaraff[k]][lv]; /\* Value of the modality of Tvaraff[k] *\/ */
- vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])];
+ /* vlv= nbcode[Tvaraff[k]][codtabm(k1,TnsdVar[Tvaraff[k]])]; */
kl++;
- sprintf(gplotcondition+strlen(gplotcondition),"$%d==%d && $%d==%d " ,kl,Tvaraff[k], kl+1, nbcode[Tvaraff[k]][lv]);
+ /* sprintf(gplotcondition+strlen(gplotcondition),"$%d==%d && $%d==%d " ,kl,Tvaraff[k], kl+1, nbcode[Tvaraff[k]][lv]); */
+ sprintf(gplotcondition+strlen(gplotcondition),"$%d==%d && $%d==%d " ,kl,lv, kl+1, vlv );
kl++;
- if(k <cptcoveff && cptcoveff>1)
+ if(k <cptcovs && cptcovs>1)
sprintf(gplotcondition+strlen(gplotcondition)," && ");
}
strcpy(gplotcondition+strlen(gplotcondition),")");
}else{
fprintf(ficgp,",\\\n '' ");
}
- if(cptcoveff ==0){ /* No covariate */
+ /* if(cptcoveff ==0){ /\* No covariate *\/ */
+ if(cptcovs ==0){ /* No covariate */
ioffset=2; /* Age is in 2 */
/*# yearproj age p11 p21 p31 p.1 p12 p22 p32 p.2 p13 p23 p33 p.3 p14 p24 p34 p.4*/
/*# 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 */
fprintf(ficgp,"#Number of graphics: first is logit, 2nd is probabilities, third is incidences per year\n");
fprintf(ficgp,"#model=1+age+%s \n",model);
fprintf(ficgp,"# Type of graphic ng=%d\n",ng);
- fprintf(ficgp,"# k1=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */
+ /* fprintf(ficgp,"# k1=1 to 2^%d=%d\n",cptcoveff,m);/\* to be checked *\/ */
+ fprintf(ficgp,"# k1=1 to 2^%d=%d\n",cptcovs,m);/* to be checked */
/* for(k1=1; k1 <=m; k1++) /\* For each combination of covariate *\/ */
for(nres=1; nres <= nresult; nres++){ /* For each resultline */
/* k1=nres; */
/* date2dmy(dateintmean,&jintmean,&mintmean,&aintmean); */
/* date2dmy(dateprojd,&jprojd, &mprojd, &anprojd); */
/* date2dmy(dateprojf,&jprojf, &mprojf, &anprojf); */
- i1=pow(2,cptcoveff);
- if (cptcovn < 1){i1=1;}
+ /* i1=pow(2,cptcoveff); */
+ /* if (cptcovn < 1){i1=1;} */
fprintf(ficresf,"# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jintmean,mintmean,aintmean,dateintmean,dateprev1,dateprev2);
fprintf(ficresf,"#****** Routine prevforecast **\n");
/* if (h==(int)(YEARM*yearp)){ */
- for(nres=1; nres <= nresult; nres++) /* For each resultline */
- for(k=1; k<=i1;k++){ /* We want to find the combination k corresponding to the values of the dummies given in this resut line (to be cleaned one day) */
- if(i1 != 1 && TKresult[nres]!= k)
- continue;
- if(invalidvarcomb[k]){
- printf("\nCombination (%d) projection ignored because no cases \n",k);
- continue;
- }
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ k=TKresult[nres];
+ if(TKresult[nres]==0) k=1; /* To be checked for noresult */
+ /* for(k=1; k<=i1;k++){ /\* We want to find the combination k corresponding to the values of the dummies given in this resut line (to be cleaned one day) *\/ */
+ /* if(i1 != 1 && TKresult[nres]!= k) */
+ /* continue; */
+ /* if(invalidvarcomb[k]){ */
+ /* printf("\nCombination (%d) projection ignored because no cases \n",k); */
+ /* continue; */
+ /* } */
fprintf(ficresf,"\n#****** hpijx=probability over h years, hp.jx is weighted by observed prev \n#");
- for(j=1;j<=cptcoveff;j++) {
- /* fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,Tvaraff[j])]); */
- fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- fprintf(ficresf," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ for(j=1;j<=cptcovs;j++){
+ /* for(j=1;j<=cptcoveff;j++) { */
+ /* /\* fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,Tvaraff[j])]); *\/ */
+ /* fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */
+ /* } */
+ /* for (k4=1; k4<= nsq; k4++){ /\* For each selected (single) quantitative value *\/ */
+ /* fprintf(ficresf," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); */
+ /* } */
+ fprintf(ficresf," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]);
}
+
fprintf(ficresf," yearproj age");
for(j=1; j<=nlstate+ndeath;j++){
for(i=1; i<=nlstate;i++)
}
}
fprintf(ficresf,"\n");
- for(j=1;j<=cptcoveff;j++)
+ /* for(j=1;j<=cptcoveff;j++) */
+ for(j=1;j<=cptcovs;j++)
+ fprintf(ficresf,"%d %lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]);
/* fprintf(ficresf,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,Tvaraff[j])]); /\* Tvaraff not correct *\/ */
- fprintf(ficresf,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); /* TnsdVar[Tvaraff] correct */
+ /* fprintf(ficresf,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); /\* TnsdVar[Tvaraff] correct *\/ */
fprintf(ficresf,"%.f %.f ",anprojd+yearp,agec+h*hstepm/YEARM*stepm);
for(j=1; j<=nlstate+ndeath;j++) {
/* if(jintmean==0) jintmean=1; */
/* if(mintmean==0) jintmean=1; */
- i1=pow(2,cptcoveff);
- if (cptcovn < 1){i1=1;}
+ /* i1=pow(2,cptcoveff); */
+ /* if (cptcovn < 1){i1=1;} */
fprintf(ficresfb,"# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jintmean,mintmean,aintmean,dateintmean,dateprev1,dateprev2);
printf("# Mean day of interviews %.lf/%.lf/%.lf (%.2f) between %.2f and %.2f \n",jintmean,mintmean,aintmean,dateintmean,dateprev1,dateprev2);
fprintf(ficresfb,"#****** Routine prevbackforecast **\n");
- for(nres=1; nres <= nresult; nres++) /* For each resultline */
- for(k=1; k<=i1;k++){
- if(i1 != 1 && TKresult[nres]!= k)
- continue;
- if(invalidvarcomb[k]){
- printf("\nCombination (%d) projection ignored because no cases \n",k);
- continue;
- }
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ k=TKresult[nres];
+ if(TKresult[nres]==0) k=1; /* To be checked for noresult */
+ /* for(k=1; k<=i1;k++){ */
+ /* if(i1 != 1 && TKresult[nres]!= k) */
+ /* continue; */
+ /* if(invalidvarcomb[k]){ */
+ /* printf("\nCombination (%d) projection ignored because no cases \n",k); */
+ /* continue; */
+ /* } */
fprintf(ficresfb,"\n#****** hbijx=probability over h years, hb.jx is weighted by observed prev \n#");
- for(j=1;j<=cptcoveff;j++) {
- fprintf(ficresfb," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]);
- }
- for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
- fprintf(ficresf," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
+ for(j=1;j<=cptcovs;j++){
+ /* for(j=1;j<=cptcoveff;j++) { */
+ /* fprintf(ficresfb," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */
+ /* } */
+ fprintf(ficresfb," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]);
}
+ /* fprintf(ficrespij,"******\n"); */
+ /* for (k4=1; k4<= nsq; k4++){ /\* For each selected (single) quantitative value *\/ */
+ /* fprintf(ficresfb," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]); */
+ /* } */
fprintf(ficresfb," yearbproj age");
for(j=1; j<=nlstate+ndeath;j++){
for(i=1; i<=nlstate;i++)
}
}
fprintf(ficresfb,"\n");
- for(j=1;j<=cptcoveff;j++)
- fprintf(ficresfb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]);
+ /* for(j=1;j<=cptcoveff;j++) */
+ for(j=1;j<=cptcovs;j++)
+ fprintf(ficresfb,"%d %lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]);
+ /* fprintf(ficresfb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */
fprintf(ficresfb,"%.f %.f ",anbackd+yearp,agec-h*hstepm/YEARM*stepm);
for(i=1; i<=nlstate+ndeath;i++) {
ppij=0.;ppi=0.;
if (strlen(modelsav) >1){ /* V2 +V3 +V4 +V6 +V7 +V6*V2 +V7*V2 +V6*V3 +V7*V3 +V6*V4 +V7*V4 +age*V2 +age*V3 +age*V4 +age*V6 +age*V7 +age*V6*V2 +V7*V2 +age*V6*V3 +age*V7*V3 +age*V6*V4 +age*V7*V4 */
j=nbocc(modelsav,'+'); /**< j=Number of '+' */
j1=nbocc(modelsav,'*'); /**< j1=Number of '*' */
- cptcovs=j+1-j1; /**< Number of simple covariates V1 +V1*age +V3 +V3*V4 +age*age => V1 + V3 =4+1-3=2 */
+ cptcovs=0; /**< Number of simple covariates V1 +V1*age +V3 +V3*V4 +age*age => V1 + V3 =4+1-3=2 Wrong */
cptcovt= j+1; /* Number of total covariates in the model, not including
* cst, age and age*age
* V1+V1*age+ V3 + V3*V4+age*age=> 3+1=4*/
Tvar[k]=0; Tprod[k]=0; Tposprod[k]=0;
}
cptcovage=0;
+
+ /* First loop in order to calculate */
+ /* for age*VN*Vm
+ * Provides, Typevar[k], Tage[cptcovage], existcomb[n][m], FixedV[ncovcolt+k12]
+ * Tprod[k1]=k Tposprod[k]=k1; Tvard[k1][1] =m;
+ */
+ /* Needs FixedV[Tvardk[k][1]] */
+ /* For others:
+ * Sets Typevar[k];
+ * Tvar[k]=ncovcol+nqv+ntv+nqtv+k11;
+ * Tposprod[k]=k11;
+ * Tprod[k11]=k;
+ * Tvardk[k][1] =m;
+ * Needs FixedV[Tvardk[k][1]] == 0
+ */
+
for(k=1; k<=cptcovt;k++){ /* Loop on total covariates of the model line */
cutl(stra,strb,modelsav,'+'); /* keeps in strb after the first '+' cutl from left to right
modelsav==V2+V1+V5*age+V4+V3*age strb=V3*age stra=V2+V1V5*age+V4 */ /* <model> "V5+V4+V3+V4*V3+V5*age+V1*age+V1" strb="V5" stra="V4+V3+V4*V3+V5*age+V1*age+V1" */
strcat(strb,stre);
strcpy(strd,strb); /* in order for strd to not be "age" for next test (will be Vn*Vm */
}
- printf("DEBUG FIXED k=%d, Tage[k]=%d, Tvar[Tage[k]=%d,FixedV[Tvar[Tage[k]]]=%d\n",k,Tage[k],Tvar[Tage[k]],FixedV[Tvar[Tage[k]]]);
- FixedV[Tvar[Tage[k]]]=0; /* HERY not sure */
+ /* printf("DEBUG FIXED k=%d, Tage[k]=%d, Tvar[Tage[k]=%d,FixedV[Tvar[Tage[k]]]=%d\n",k,Tage[k],Tvar[Tage[k]],FixedV[Tvar[Tage[k]]]); */
+ /* FixedV[Tvar[Tage[k]]]=0; /\* HERY not sure if V7*V4*age Fixed might not exist yet*\/ */
}else{ /* strc=Vn*Vm (and strd=age) and should be strb=Vn*Vm but want to keep original strb double product */
strcpy(stre,strb); /* save full b in stre */
strcpy(strb,strc); /* save short c in new short b for next block strb=Vn*Vm*/
Tvardk[k][1] =m; /* m 1 for V1*/
Tvard[k1][2] =n; /* n 4 for V4*/
Tvardk[k][2] =n; /* n 4 for V4*/
-/* Tvar[Tage[cptcovage]]=k1;*/ /* Tvar[6=age*V3*V2]=9 (new fixed covariate) */
+/* Tvar[Tage[cptcovage]]=k1;*/ /* Tvar[6=age*V3*V2]=9 (new fixed covariate) */ /* We don't know about Fixed yet HERE */
if( FixedV[Tvardk[k][1]] == 0 && FixedV[Tvardk[k][2]] == 0){ /* If the product is a fixed covariate then we feed the new column with Vn*Vm */
for (i=1; i<=lastobs;i++){/* For fixed product */
/* Computes the new covariate which is a product of
/*printf("a=%s b=%s sav=%s\n", stra,strb,modelsav);
scanf("%d",i);*/
} /* end of loop + on total covariates */
+
+
} /* end if strlen(modelsave == 0) age*age might exist */
} /* end if strlen(model == 0) */
cptcovs=cptcovt - cptcovdageprod - cptcovprod;/**< Number of simple covariates V1 +V1*age +V3 +V3*V4 +age*age + age*v4*V3=> V1 + V3 =4+1-3=2 */
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<=NCOVMAX; k++){ Fixed[k]=0; Dummy[k]=0;}
for(k=1;k<=NCOVMAX; k++){TvarFind[k]=0; TvarVind[k]=0;}
+
+
+ /* Second loop for calculating Fixed[k], Dummy[k]*/
+
+
for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0,ncovva=0,ncovvta=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0, ncovvt=0;k<=cptcovt; k++){ /* or cptocvt loop on k from model */
if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */
Fixed[k]= 0;
/* double ***mobaverage; */
double wald;
- char line[MAXLINE];
+ char line[MAXLINE], linetmp[MAXLINE];
char path[MAXLINE],pathc[MAXLINE],pathcd[MAXLINE],pathtot[MAXLINE];
char modeltemp[MAXLINE];
}else
break;
}
- if((num_filled=sscanf(line,"model=1+age%[^.\n]", model)) !=EOF){
+ if((num_filled=sscanf(line,"model=%[^.\n]", model)) !=EOF){ /* Every character after model but dot and return */
+ if (num_filled != 1){
+ printf("ERROR %d: Model should be at minimum 'model=1+age+' instead of '%s'\n",num_filled, line);
+ fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age+' instead of '%s'\n",num_filled, line);
+ model[0]='\0';
+ goto end;
+ }else{
+ trimbtab(linetmp,line); /* Trims multiple blanks in line */
+ strcpy(line, linetmp);
+ }
+ }
+ if((num_filled=sscanf(line,"model=1+age%[^.\n]", model)) !=EOF){ /* Every character after 1+age but dot and return */
if (num_filled != 1){
printf("ERROR %d: Model should be at minimum 'model=1+age+' instead of '%s'\n",num_filled, line);
fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age+' instead of '%s'\n",num_filled, line);
Tvard=imatrix(1,NCOVMAX,1,2); /* n=Tvard[k1][1] and m=Tvard[k1][2] gives the couple n,m of the k1 th product Vn*Vm
* For V3*V2 (in V2+V1+V1*V4+age*V3+V3*V2), V3*V2 position is 2nd.
* Tvard[k1=2][1]=3 (V3) Tvard[k1=2][2]=2(V2) */
- Tvardk=imatrix(-1,NCOVMAX,1,2);
+ Tvardk=imatrix(0,NCOVMAX,1,2);
Tage=ivector(1,NCOVMAX); /* Gives the covariate id of covariates associated with age: V2 + V1 + age*V4 + V3*age
4 covariates (3 plus signs)
Tage[1=V3*age]= 4; Tage[2=age*V4] = 3
pstamp(ficreseij);
- i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */
- if (cptcovn < 1){i1=1;}
+ /* i1=pow(2,cptcoveff); /\* Number of combination of dummy covariates *\/ */
+ /* if (cptcovn < 1){i1=1;} */
- for(nres=1; nres <= nresult; nres++) /* For each resultline */
- for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
- if(i1 != 1 && TKresult[nres]!= k)
- continue;
+ for(nres=1; nres <= nresult; nres++){ /* For each resultline */
+ /* for(k=1; k<=i1;k++){ /\* For any combination of dummy covariates, fixed and varying *\/ */
+ /* if(i1 != 1 && TKresult[nres]!= k) */
+ /* continue; */
fprintf(ficreseij,"\n#****** ");
printf("\n#****** ");
- for(j=1;j<=cptcoveff;j++) {
- fprintf(ficreseij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]);
- printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]);
+ for(j=1;j<=cptcovs;j++){
+ /* for(j=1;j<=cptcoveff;j++) { */
+ /* fprintf(ficreseij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */
+ fprintf(ficreseij," V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]);
+ printf(" V%d=%lg ",Tvresult[nres][j],TinvDoQresult[nres][Tvresult[nres][j]]);
+ /* printf("V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,TnsdVar[Tvaraff[j]])]); */
}
for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
printf(" V%d=%lg ",TvarsQ[j], TinvDoQresult[nres][TvarsQ[j]]); /* TvarsQ[j] gives the name of the jth quantitative (fixed or time v) */
free_vector(weight,firstobs,lastobs);
- free_imatrix(Tvardk,-1,NCOVMAX,1,2);
+ free_imatrix(Tvardk,0,NCOVMAX,1,2);
free_imatrix(Tvard,1,NCOVMAX,1,2);
free_imatrix(s,1,maxwav+1,firstobs,lastobs);
free_matrix(anint,1,maxwav,firstobs,lastobs);
git://henry.ined.fr / .git / commitdiff