SWZ: replace tabs with spaces
parent
8b782e7b00
commit
60a3c2cad4
|
@ -562,7 +562,7 @@ int ReadBaseTransitionMatrices_SV(FILE *f_in, TMarkovStateVariable* sv, char *he
|
|||
ConvertBaseTransitionMatrix(sv->Q,Q,sv->nlags_encoded);
|
||||
|
||||
/* Free Q */
|
||||
FreeMatrix(Q);
|
||||
FreeMatrix(Q);
|
||||
|
||||
/* // Update ansi-c*/
|
||||
if (!Update_B_from_Q_SV(sv))
|
||||
|
@ -728,42 +728,42 @@ int ReadBaseTransitionMatricesFlat_SV(FILE *f_in, TMarkovStateVariable *sv)
|
|||
if (sv->n_state_variables > 1)
|
||||
{
|
||||
for (i=0; i < sv->n_state_variables; i++)
|
||||
if (!ReadBaseTransitionMatricesFlat_SV(f_in,sv->state_variable[i]))
|
||||
return 0;
|
||||
if (!ReadBaseTransitionMatricesFlat_SV(f_in,sv->state_variable[i]))
|
||||
return 0;
|
||||
}
|
||||
else
|
||||
{
|
||||
// Read transition matrix
|
||||
Q=CreateMatrix(sv->nbasestates,sv->nbasestates);
|
||||
for (j=0; j < ColM(Q); j++)
|
||||
for (i=0; i < RowM(Q); i++)
|
||||
if (fscanf(f_in," %lf ",&ElementM(Q,i,j)) != 1)
|
||||
{
|
||||
FreeMatrix(Q);
|
||||
return 0;
|
||||
}
|
||||
for (i=0; i < RowM(Q); i++)
|
||||
if (fscanf(f_in," %lf ",&ElementM(Q,i,j)) != 1)
|
||||
{
|
||||
FreeMatrix(Q);
|
||||
return 0;
|
||||
}
|
||||
|
||||
// Scale the columns of Q - loose requirement on sumation to one
|
||||
for (j=sv->nbasestates-1; j >= 0; j--)
|
||||
{
|
||||
for (sum=0.0, i=sv->nbasestates-1; i >= 0; i--)
|
||||
if (ElementM(Q,i,j) < 0.0)
|
||||
{
|
||||
FreeMatrix(Q);
|
||||
dw_UserError("Transition matrix can not have negative elements.");
|
||||
return 0;
|
||||
}
|
||||
else
|
||||
sum+=ElementM(Q,i,j);
|
||||
if (fabs(sum-1.0) > 1.0e-4)
|
||||
{
|
||||
FreeMatrix(Q);
|
||||
dw_UserError("Transition matrix columns must sum to one.");
|
||||
return 0;
|
||||
}
|
||||
for (sum=1.0/sum, i=sv->nbasestates-1; i >= 0; i--)
|
||||
ElementM(Q,i,j)*=sum;
|
||||
}
|
||||
{
|
||||
for (sum=0.0, i=sv->nbasestates-1; i >= 0; i--)
|
||||
if (ElementM(Q,i,j) < 0.0)
|
||||
{
|
||||
FreeMatrix(Q);
|
||||
dw_UserError("Transition matrix can not have negative elements.");
|
||||
return 0;
|
||||
}
|
||||
else
|
||||
sum+=ElementM(Q,i,j);
|
||||
if (fabs(sum-1.0) > 1.0e-4)
|
||||
{
|
||||
FreeMatrix(Q);
|
||||
dw_UserError("Transition matrix columns must sum to one.");
|
||||
return 0;
|
||||
}
|
||||
for (sum=1.0/sum, i=sv->nbasestates-1; i >= 0; i--)
|
||||
ElementM(Q,i,j)*=sum;
|
||||
}
|
||||
|
||||
// Convert base transition matrix to full transition matrix.
|
||||
ConvertBaseTransitionMatrix(sv->Q,Q,sv->nlags_encoded);
|
||||
|
@ -773,10 +773,10 @@ int ReadBaseTransitionMatricesFlat_SV(FILE *f_in, TMarkovStateVariable *sv)
|
|||
|
||||
// Update
|
||||
if (!Update_B_from_Q_SV(sv))
|
||||
{
|
||||
dw_UserError("Transition matrices do not satisfy restrictions");
|
||||
return 0;
|
||||
}
|
||||
{
|
||||
dw_UserError("Transition matrices do not satisfy restrictions");
|
||||
return 0;
|
||||
}
|
||||
|
||||
return 1;
|
||||
}
|
||||
|
|
|
@ -2433,13 +2433,13 @@ TMatrix forecast_base(TMatrix forecast, int horizon, TVector initial, TVector *s
|
|||
{
|
||||
ProductVM(y,x,Aplus[S[t]]);
|
||||
if (shocks)
|
||||
{
|
||||
for (i=p->nvars-1; i >= 0; i--)
|
||||
ElementV(y,i)+=ElementV(shocks[t],i)/sqrt(p->Zeta[i][p->var_states[i][S[t]]]);
|
||||
}
|
||||
{
|
||||
for (i=p->nvars-1; i >= 0; i--)
|
||||
ElementV(y,i)+=ElementV(shocks[t],i)/sqrt(p->Zeta[i][p->var_states[i][S[t]]]);
|
||||
}
|
||||
ProductInverseVM(y,y,A0[S[t]]);
|
||||
for (i=p->nvars-1; i >= 0; i--)
|
||||
ElementM(forecast,t,i)=ElementV(y,i);
|
||||
ElementM(forecast,t,i)=ElementV(y,i);
|
||||
|
||||
memmove(pElementV(x)+p->nvars,pElementV(x),(p->nlags-1)*p->nvars*sizeof(PRECISION));
|
||||
memcpy(pElementV(x),pElementV(y),p->nvars*sizeof(PRECISION));
|
||||
|
@ -2487,23 +2487,23 @@ TMatrix forecast_base(TMatrix forecast, int horizon, TVector initial, TVector *s
|
|||
/* Pxi1[k]=CreateVector(statespace->nstates); */
|
||||
/* SPxi[k]=CreateVector(statespace->nstates); */
|
||||
/* for (i=statespace->zeta_modulus-1; i >= 0; i--) */
|
||||
/* { */
|
||||
/* IEz[k][i]=CreateVector(statespace->nz); */
|
||||
/* IEzz[k][i]=CreateMatrix(statespace->nz,statespace->nz); */
|
||||
/* ISEz[k][i]=CreateVector(statespace->nz); */
|
||||
/* } */
|
||||
/* { */
|
||||
/* IEz[k][i]=CreateVector(statespace->nz); */
|
||||
/* IEzz[k][i]=CreateMatrix(statespace->nz,statespace->nz); */
|
||||
/* ISEz[k][i]=CreateVector(statespace->nz); */
|
||||
/* } */
|
||||
/* for (i=statespace->nstates-1; i >= 0; i--) */
|
||||
/* { */
|
||||
/* Ez1[k][i]=CreateVector(statespace->nz); */
|
||||
/* Ezz1[k][i]=CreateMatrix(statespace->nz,statespace->nz); */
|
||||
/* Ez[k][i]=CreateVector(statespace->nz); */
|
||||
/* Ezz[k][i]=CreateMatrix(statespace->nz,statespace->nz); */
|
||||
/* SEz[k][i]=CreateVector(statespace->nz); */
|
||||
/* } */
|
||||
/* { */
|
||||
/* Ez1[k][i]=CreateVector(statespace->nz); */
|
||||
/* Ezz1[k][i]=CreateMatrix(statespace->nz,statespace->nz); */
|
||||
/* Ez[k][i]=CreateVector(statespace->nz); */
|
||||
/* Ezz[k][i]=CreateMatrix(statespace->nz,statespace->nz); */
|
||||
/* SEz[k][i]=CreateVector(statespace->nz); */
|
||||
/* } */
|
||||
/* } */
|
||||
|
||||
/* ConditionalFilter(0,h-1,y,statespace->Ez[t0],statespace->Ezz[t0],model->V[t0],model->sv->Q, */
|
||||
/* Pxi,Pxi1,IEz,IEzz,Ez1,Ezz1,Ez,Ezz,statespace); */
|
||||
/* Pxi,Pxi1,IEz,IEzz,Ez1,Ezz1,Ez,Ezz,statespace); */
|
||||
|
||||
/* SmoothProbabilities_MSStateSpace(0,h-1,SPxi,Pxi,Pxi1,model->sv->Q); */
|
||||
|
||||
|
@ -2519,34 +2519,34 @@ TMatrix forecast_base(TMatrix forecast, int horizon, TVector initial, TVector *s
|
|||
/* { */
|
||||
/* F=dw_CreateArray_vector(h); */
|
||||
/* for (k=h-1; k >= 0; k--) */
|
||||
/* F[k]=CreateVector(statespace->ny); */
|
||||
/* F[k]=CreateVector(statespace->ny); */
|
||||
/* } */
|
||||
|
||||
/* for (k=h-1; k >= 0; k--) */
|
||||
/* { */
|
||||
/* InitializeVector(F[k],0.0); */
|
||||
/* if (statespace->zeta_modulus > statespace->nbasestates) */
|
||||
/* { */
|
||||
/* IntegrateStatesSingle(SPzeta,SPxi[k],statespace->zeta_modulus,statespace->nbasestates,2); */
|
||||
/* IntegrateStatesSingleV(z,SPzeta,ISEz[k],statespace->nbasestates,statespace->zeta_modulus/statespace->nbasestates,2); */
|
||||
/* IntegrateStatesSingle(SPs,SPzeta,statespace->nbasestates,statespace->zeta_modulus/statespace->nbasestates,2); */
|
||||
/* for (s=statespace->nbasestates-1; s >= 0; s--) */
|
||||
/* { */
|
||||
/* ProductMV(u,statespace->H[s],z[s]); */
|
||||
/* AddVV(u,statespace->a[s],u); */
|
||||
/* LinearCombinationV(F[k],1.0,F[k],ElementV(SPs,s),u); */
|
||||
/* } */
|
||||
/* } */
|
||||
/* { */
|
||||
/* IntegrateStatesSingle(SPzeta,SPxi[k],statespace->zeta_modulus,statespace->nbasestates,2); */
|
||||
/* IntegrateStatesSingleV(z,SPzeta,ISEz[k],statespace->nbasestates,statespace->zeta_modulus/statespace->nbasestates,2); */
|
||||
/* IntegrateStatesSingle(SPs,SPzeta,statespace->nbasestates,statespace->zeta_modulus/statespace->nbasestates,2); */
|
||||
/* for (s=statespace->nbasestates-1; s >= 0; s--) */
|
||||
/* { */
|
||||
/* ProductMV(u,statespace->H[s],z[s]); */
|
||||
/* AddVV(u,statespace->a[s],u); */
|
||||
/* LinearCombinationV(F[k],1.0,F[k],ElementV(SPs,s),u); */
|
||||
/* } */
|
||||
/* } */
|
||||
/* else */
|
||||
/* { */
|
||||
/* IntegrateStatesSingle(SPs,SPxi[k],statespace->nbasestates,statespace->zeta_modulus,2); */
|
||||
/* for (s=statespace->nbasestates-1; s >= 0; s--) */
|
||||
/* { */
|
||||
/* ProductMV(u,statespace->H[s],ISEz[k][s]); */
|
||||
/* AddVV(u,statespace->a[s],u); */
|
||||
/* LinearCombinationV(F[k],1.0,F[k],ElementV(SPs,s),u); */
|
||||
/* } */
|
||||
/* } */
|
||||
/* { */
|
||||
/* IntegrateStatesSingle(SPs,SPxi[k],statespace->nbasestates,statespace->zeta_modulus,2); */
|
||||
/* for (s=statespace->nbasestates-1; s >= 0; s--) */
|
||||
/* { */
|
||||
/* ProductMV(u,statespace->H[s],ISEz[k][s]); */
|
||||
/* AddVV(u,statespace->a[s],u); */
|
||||
/* LinearCombinationV(F[k],1.0,F[k],ElementV(SPs,s),u); */
|
||||
/* } */
|
||||
/* } */
|
||||
/* } */
|
||||
|
||||
/* // Clean up */
|
||||
|
@ -3590,19 +3590,19 @@ TMatrix* MakeA0_All(TMatrix *A0, T_VAR_Parameters *p)
|
|||
if (!A0)
|
||||
{
|
||||
if (!(A0=dw_CreateArray_matrix(p->nstates)))
|
||||
return (TMatrix*)NULL;
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
else
|
||||
if ((dw_DimA(A0) != p->nstates))
|
||||
{
|
||||
dw_Error(SIZE_ERR);
|
||||
return (TMatrix*)NULL;
|
||||
dw_Error(SIZE_ERR);
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
for (s=p->nstates-1; s >= 0; s--)
|
||||
if (!(A0[s]=MakeA0(A0[s],s,p)))
|
||||
{
|
||||
if (A0_in != A0) dw_FreeArray(A0);
|
||||
return (TMatrix*)NULL;
|
||||
if (A0_in != A0) dw_FreeArray(A0);
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
return A0;
|
||||
}
|
||||
|
@ -3644,19 +3644,19 @@ TMatrix* MakeAplus_All(TMatrix *Aplus, T_VAR_Parameters *p)
|
|||
if (!Aplus)
|
||||
{
|
||||
if (!(Aplus=dw_CreateArray_matrix(p->nstates)))
|
||||
return (TMatrix*)NULL;
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
else
|
||||
if ((dw_DimA(Aplus) != p->nstates))
|
||||
{
|
||||
dw_Error(SIZE_ERR);
|
||||
return (TMatrix*)NULL;
|
||||
dw_Error(SIZE_ERR);
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
for (s=p->nstates-1; s >= 0; s--)
|
||||
if (!(Aplus[s]=MakeAplus(Aplus[s],s,p)))
|
||||
{
|
||||
if (Aplus_in != Aplus) dw_FreeArray(Aplus);
|
||||
return (TMatrix*)NULL;
|
||||
if (Aplus_in != Aplus) dw_FreeArray(Aplus);
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
return Aplus;
|
||||
}
|
||||
|
@ -3694,19 +3694,19 @@ TMatrix* MakeZeta_All(TMatrix *Zeta, T_VAR_Parameters *p)
|
|||
if (!Zeta)
|
||||
{
|
||||
if (!(Zeta=dw_CreateArray_matrix(p->nstates)))
|
||||
return (TMatrix*)NULL;
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
else
|
||||
if ((dw_DimA(Zeta) != p->nstates))
|
||||
{
|
||||
dw_Error(SIZE_ERR);
|
||||
return (TMatrix*)NULL;
|
||||
dw_Error(SIZE_ERR);
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
for (s=p->nstates-1; s >= 0; s--)
|
||||
if (!(Zeta[s]=MakeZeta(Zeta[s],s,p)))
|
||||
{
|
||||
if (Zeta_in != Zeta) dw_FreeArray(Zeta);
|
||||
return (TMatrix*)NULL;
|
||||
if (Zeta_in != Zeta) dw_FreeArray(Zeta);
|
||||
return (TMatrix*)NULL;
|
||||
}
|
||||
return Zeta;
|
||||
}
|
||||
|
|
|
@ -580,36 +580,36 @@ int Read_VAR_ParametersFlat(FILE *f_in, TStateModel *model)
|
|||
{
|
||||
A0[s]=CreateMatrix(p->nvars,p->nvars);
|
||||
for (j=0; j < p->nvars; j++)
|
||||
for (i=0; i < p->nvars; i++)
|
||||
if (fscanf(f_in," %lf ",&ElementM(A0[s],i,j)) != 1)
|
||||
goto ERROR;
|
||||
for (i=0; i < p->nvars; i++)
|
||||
if (fscanf(f_in," %lf ",&ElementM(A0[s],i,j)) != 1)
|
||||
goto ERROR;
|
||||
|
||||
Aplus[s]=CreateMatrix(p->npre,p->nvars);
|
||||
for (j=0; j < p->nvars; j++)
|
||||
for (i=0; i < p->npre; i++)
|
||||
if (fscanf(f_in," %lf ",&ElementM(Aplus[s],i,j)) != 1)
|
||||
goto ERROR;
|
||||
for (i=0; i < p->npre; i++)
|
||||
if (fscanf(f_in," %lf ",&ElementM(Aplus[s],i,j)) != 1)
|
||||
goto ERROR;
|
||||
|
||||
Zeta[s]=CreateVector(p->nvars);
|
||||
for (j=0; j < p->nvars; j++)
|
||||
if (fscanf(f_in," %lf ",&ElementV(Zeta[s],j)) != 1)
|
||||
goto ERROR;
|
||||
else
|
||||
if (ElementV(Zeta[s],j) < 0.0)
|
||||
goto ERROR;
|
||||
if (fscanf(f_in," %lf ",&ElementV(Zeta[s],j)) != 1)
|
||||
goto ERROR;
|
||||
else
|
||||
if (ElementV(Zeta[s],j) < 0.0)
|
||||
goto ERROR;
|
||||
}
|
||||
|
||||
// Set A0, Aplus, and Zeta
|
||||
for (j=0; j < p->nvars; j++)
|
||||
for (s=0; s < p->nstates; s++)
|
||||
{
|
||||
for (i=0; i < p->nvars; i++)
|
||||
ElementV(p->A0[j][p->coef_states[j][s]],i)=ElementM(A0[s],i,j);
|
||||
for (i=0; i < p->nvars; i++)
|
||||
ElementV(p->A0[j][p->coef_states[j][s]],i)=ElementM(A0[s],i,j);
|
||||
|
||||
for (i=0; i < p->npre; i++)
|
||||
ElementV(p->Aplus[j][p->coef_states[j][s]],i)=ElementM(Aplus[s],i,j);
|
||||
for (i=0; i < p->npre; i++)
|
||||
ElementV(p->Aplus[j][p->coef_states[j][s]],i)=ElementM(Aplus[s],i,j);
|
||||
|
||||
p->Zeta[j][p->var_states[j][s]]=ElementV(Zeta[s],j);
|
||||
p->Zeta[j][p->var_states[j][s]]=ElementV(Zeta[s],j);
|
||||
}
|
||||
|
||||
// Update b0, bplus, lambda, psi
|
||||
|
|
|
@ -62,49 +62,49 @@ int forecast_percentile(FILE *f_out, TVector percentiles, int draws, FILE *poste
|
|||
{
|
||||
// Read parameters and push them into model
|
||||
if (!posterior_file)
|
||||
done=1;
|
||||
done=1;
|
||||
else
|
||||
if (!ReadBaseTransitionMatricesFlat(posterior_file,model) || !Read_VAR_ParametersFlat(posterior_file,model))
|
||||
{
|
||||
done=2;
|
||||
printf("total posterior draws processed - %d\n",i);
|
||||
}
|
||||
else
|
||||
if (i++ == n)
|
||||
{
|
||||
printf("%d posterior draws processed\n",i);
|
||||
n+=1000;
|
||||
}
|
||||
if (!ReadBaseTransitionMatricesFlat(posterior_file,model) || !Read_VAR_ParametersFlat(posterior_file,model))
|
||||
{
|
||||
done=2;
|
||||
printf("total posterior draws processed - %d\n",i);
|
||||
}
|
||||
else
|
||||
if (i++ == n)
|
||||
{
|
||||
printf("%d posterior draws processed\n",i);
|
||||
n+=1000;
|
||||
}
|
||||
|
||||
if (done != 2)
|
||||
{
|
||||
// Get filtered probability at time T
|
||||
for (j=p->nstates-1; j >= 0; j--)
|
||||
ElementV(init_prob,j)=ProbabilityStateConditionalCurrent(j,T,model);
|
||||
if (done != 2)
|
||||
{
|
||||
// Get filtered probability at time T
|
||||
for (j=p->nstates-1; j >= 0; j--)
|
||||
ElementV(init_prob,j)=ProbabilityStateConditionalCurrent(j,T,model);
|
||||
|
||||
for (k=draws; k > 0; k--)
|
||||
{
|
||||
// Draw time T regime
|
||||
m=DrawDiscrete(init_prob);
|
||||
for (k=draws; k > 0; k--)
|
||||
{
|
||||
// Draw time T regime
|
||||
m=DrawDiscrete(init_prob);
|
||||
|
||||
// Draw regimes from time T+1 through T+h inclusive
|
||||
for (j=0; j < h; j++)
|
||||
{
|
||||
ColumnVector(prob,model->sv->Q,m);
|
||||
S[j]=m=DrawDiscrete(prob);
|
||||
}
|
||||
// Draw regimes from time T+1 through T+h inclusive
|
||||
for (j=0; j < h; j++)
|
||||
{
|
||||
ColumnVector(prob,model->sv->Q,m);
|
||||
S[j]=m=DrawDiscrete(prob);
|
||||
}
|
||||
|
||||
// Draw shocks
|
||||
for (j=h-1; j >= 0; j--) dw_NormalVector(shocks[j]); // InitializeVector(shocks[i],0.0);
|
||||
// Draw shocks
|
||||
for (j=h-1; j >= 0; j--) dw_NormalVector(shocks[j]); // InitializeVector(shocks[i],0.0);
|
||||
|
||||
// Compute forecast
|
||||
if (!forecast_base(forecast,h,initial,shocks,S,model))
|
||||
goto ERROR_EXIT;
|
||||
// Compute forecast
|
||||
if (!forecast_base(forecast,h,initial,shocks,S,model))
|
||||
goto ERROR_EXIT;
|
||||
|
||||
// Accumulate impulse response
|
||||
AddMatrixObservation(forecast,histogram);
|
||||
}
|
||||
}
|
||||
// Accumulate impulse response
|
||||
AddMatrixObservation(forecast,histogram);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
for (i=0; i < DimV(percentiles); i++)
|
||||
|
@ -152,7 +152,7 @@ ERROR_EXIT:
|
|||
|
||||
*/
|
||||
int forecast_percentile_regime(FILE *f_out, TVector percentiles, int draws,
|
||||
FILE *posterior_file, int s, int T, int h, TStateModel *model)
|
||||
FILE *posterior_file, int s, int T, int h, TStateModel *model)
|
||||
{
|
||||
T_VAR_Parameters *p;
|
||||
int done=0, rtrn=0, *S, i=0, j, k, m, n=1000;
|
||||
|
@ -186,35 +186,35 @@ int forecast_percentile_regime(FILE *f_out, TVector percentiles, int draws,
|
|||
{
|
||||
// Read parameters and push them into model
|
||||
if (!posterior_file)
|
||||
done=1;
|
||||
done=1;
|
||||
else
|
||||
if (!ReadBaseTransitionMatricesFlat(posterior_file,model) || !Read_VAR_ParametersFlat(posterior_file,model))
|
||||
{
|
||||
done=2;
|
||||
printf("total posterior draws processed - %d\n",i);
|
||||
}
|
||||
else
|
||||
if (i++ == n)
|
||||
{
|
||||
printf("%d posterior draws processed\n",i);
|
||||
n+=1000;
|
||||
}
|
||||
if (!ReadBaseTransitionMatricesFlat(posterior_file,model) || !Read_VAR_ParametersFlat(posterior_file,model))
|
||||
{
|
||||
done=2;
|
||||
printf("total posterior draws processed - %d\n",i);
|
||||
}
|
||||
else
|
||||
if (i++ == n)
|
||||
{
|
||||
printf("%d posterior draws processed\n",i);
|
||||
n+=1000;
|
||||
}
|
||||
|
||||
if (done != 2)
|
||||
{
|
||||
for (k=draws; k > 0; k--)
|
||||
{
|
||||
// Draw shocks
|
||||
for (j=h-1; j >= 0; j--) dw_NormalVector(shocks[j]); // InitializeVector(shocks[i],0.0);
|
||||
{
|
||||
for (k=draws; k > 0; k--)
|
||||
{
|
||||
// Draw shocks
|
||||
for (j=h-1; j >= 0; j--) dw_NormalVector(shocks[j]); // InitializeVector(shocks[i],0.0);
|
||||
|
||||
// Compute forecast
|
||||
if (!forecast_base(forecast,h,initial,shocks,S,model))
|
||||
goto ERROR_EXIT;
|
||||
// Compute forecast
|
||||
if (!forecast_base(forecast,h,initial,shocks,S,model))
|
||||
goto ERROR_EXIT;
|
||||
|
||||
// Accumulate impulse response
|
||||
AddMatrixObservation(forecast,histogram);
|
||||
}
|
||||
}
|
||||
// Accumulate impulse response
|
||||
AddMatrixObservation(forecast,histogram);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
for (i=0; i < DimV(percentiles); i++)
|
||||
|
@ -325,11 +325,11 @@ int main(int nargs, char **args)
|
|||
|
||||
// specification filename
|
||||
if (!spec)
|
||||
sprintf(spec=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
sprintf(spec=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
|
||||
// parameter filename
|
||||
if (!parm)
|
||||
sprintf(parm=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
sprintf(parm=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
}
|
||||
|
||||
// horizon
|
||||
|
@ -340,11 +340,11 @@ int main(int nargs, char **args)
|
|||
swz_fprintf_err("No specification filename given\n");
|
||||
swz_fprintf_err("Command line syntax:\n"
|
||||
" -ft : file tag\n"
|
||||
" -fs : specification filename (est_final_<tag>.dat)\n"
|
||||
" -fp : parameters filename (specification filename)\n"
|
||||
" -fs : specification filename (est_final_<tag>.dat)\n"
|
||||
" -fp : parameters filename (specification filename)\n"
|
||||
" -fh : parameter header (Posterior mode: )\n"
|
||||
" -horizon : horizon for the forecast (12)\n"
|
||||
);
|
||||
);
|
||||
swzExit(1);
|
||||
}
|
||||
|
||||
|
@ -377,8 +377,8 @@ int main(int nargs, char **args)
|
|||
/* swzFree(out_filename); */
|
||||
/* printf("Constructing mean forecast\n"); */
|
||||
/* if (F=dw_state_space_mean_unconditional_forecast((TVector*)NULL,h,statespace->nobs,model)) */
|
||||
/* for (i=0; i < h; i++) */
|
||||
/* dw_PrintVector(f_out,F[i],"%le "); */
|
||||
/* for (i=0; i < h; i++) */
|
||||
/* dw_PrintVector(f_out,F[i],"%le "); */
|
||||
/* fclose(f_out); */
|
||||
/* return; */
|
||||
/* } */
|
||||
|
@ -390,10 +390,10 @@ int main(int nargs, char **args)
|
|||
fmt="draws_%s.dat";
|
||||
sprintf(post=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
if (!(posterior_file=fopen(post,"rt")))
|
||||
{
|
||||
printf("Unable to open draws file: %s\n",post);
|
||||
swzExit(0);
|
||||
}
|
||||
{
|
||||
printf("Unable to open draws file: %s\n",post);
|
||||
swzExit(0);
|
||||
}
|
||||
|
||||
// Get thinning factor from command line
|
||||
thin=dw_ParseInteger_String(nargs,args,"thin",1);
|
||||
|
@ -417,66 +417,66 @@ int main(int nargs, char **args)
|
|||
if ((i=dw_FindArgument_String(nargs,args,"percentiles")) == -1)
|
||||
if (dw_FindArgument_String(nargs,args,"error_bands") == -1)
|
||||
{
|
||||
percentiles=CreateVector(1);
|
||||
ElementV(percentiles,0)=0.5;
|
||||
percentiles=CreateVector(1);
|
||||
ElementV(percentiles,0)=0.5;
|
||||
}
|
||||
else
|
||||
{
|
||||
percentiles=CreateVector(3);
|
||||
ElementV(percentiles,0)=0.16; ElementV(percentiles,1)=0.5; ElementV(percentiles,2)=0.84;
|
||||
percentiles=CreateVector(3);
|
||||
ElementV(percentiles,0)=0.16; ElementV(percentiles,1)=0.5; ElementV(percentiles,2)=0.84;
|
||||
}
|
||||
else
|
||||
if ((i+1 < nargs) && dw_IsInteger(args[i+1]) && ((n=atoi(args[i+1])) > 0) && (i+1+n < nargs))
|
||||
{
|
||||
percentiles=CreateVector(n);
|
||||
for (j=0; j < n; j++)
|
||||
if (!dw_IsFloat(args[i+2+j])|| ((ElementV(percentiles,j)=atof(args[i+2+j])) <= 0.0)
|
||||
|| (ElementV(percentiles,j) >= 1.0)) break;
|
||||
if (j < n)
|
||||
{
|
||||
FreeVector(percentiles);
|
||||
printf("forecasting command line: Error parsing percentiles\n");
|
||||
swzExit(0);
|
||||
}
|
||||
percentiles=CreateVector(n);
|
||||
for (j=0; j < n; j++)
|
||||
if (!dw_IsFloat(args[i+2+j])|| ((ElementV(percentiles,j)=atof(args[i+2+j])) <= 0.0)
|
||||
|| (ElementV(percentiles,j) >= 1.0)) break;
|
||||
if (j < n)
|
||||
{
|
||||
FreeVector(percentiles);
|
||||
printf("forecasting command line: Error parsing percentiles\n");
|
||||
swzExit(0);
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
printf("forecasting command line(): Error parsing percentiles\n");
|
||||
swzExit(0);
|
||||
printf("forecasting command line(): Error parsing percentiles\n");
|
||||
swzExit(0);
|
||||
}
|
||||
|
||||
if (dw_FindArgument_String(nargs,args,"regimes") != -1)
|
||||
for (s=0; s < p->nstates; s++)
|
||||
{
|
||||
rewind(posterior_file);
|
||||
fmt="forecasts_percentiles_regime_%d_%s.prn";
|
||||
sprintf(out_filename=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 3),fmt,s,tag);
|
||||
f_out=fopen(out_filename,"wt");
|
||||
swzFree(out_filename);
|
||||
printf("Constructing percentiles for forecasts - regime %d\n",s);
|
||||
forecast_percentile_regime(f_out,percentiles,draws,posterior_file,s,p->nobs,horizon,model);
|
||||
fclose(f_out);
|
||||
rewind(posterior_file);
|
||||
fmt="forecasts_percentiles_regime_%d_%s.prn";
|
||||
sprintf(out_filename=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 3),fmt,s,tag);
|
||||
f_out=fopen(out_filename,"wt");
|
||||
swzFree(out_filename);
|
||||
printf("Constructing percentiles for forecasts - regime %d\n",s);
|
||||
forecast_percentile_regime(f_out,percentiles,draws,posterior_file,s,p->nobs,horizon,model);
|
||||
fclose(f_out);
|
||||
}
|
||||
else
|
||||
if (((s=dw_ParseInteger_String(nargs,args,"regime",-1)) >= 0) && (s < p->nstates))
|
||||
{
|
||||
fmt="forecasts_percentiles_regime_%d_%s.prn";
|
||||
sprintf(out_filename=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 3),fmt,s,tag);
|
||||
f_out=fopen(out_filename,"wt");
|
||||
swzFree(out_filename);
|
||||
printf("Constructing percentiles for forecasts - regime %d\n",s);
|
||||
forecast_percentile_regime(f_out,percentiles,draws,posterior_file,s,p->nobs,horizon,model);
|
||||
fclose(f_out);
|
||||
fmt="forecasts_percentiles_regime_%d_%s.prn";
|
||||
sprintf(out_filename=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 3),fmt,s,tag);
|
||||
f_out=fopen(out_filename,"wt");
|
||||
swzFree(out_filename);
|
||||
printf("Constructing percentiles for forecasts - regime %d\n",s);
|
||||
forecast_percentile_regime(f_out,percentiles,draws,posterior_file,s,p->nobs,horizon,model);
|
||||
fclose(f_out);
|
||||
}
|
||||
else
|
||||
{
|
||||
fmt="forecasts_percentiles_%s.prn";
|
||||
sprintf(out_filename=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
f_out=fopen(out_filename,"wt");
|
||||
swzFree(out_filename);
|
||||
printf("Constructing percentiles for forecasts - %d draws of shocks/regimes per posterior value\n",draws);
|
||||
forecast_percentile(f_out,percentiles,draws,posterior_file,p->nobs,horizon,model);
|
||||
fclose(f_out);
|
||||
fmt="forecasts_percentiles_%s.prn";
|
||||
sprintf(out_filename=(char*)swzMalloc(strlen(fmt) + strlen(tag) - 1),fmt,tag);
|
||||
f_out=fopen(out_filename,"wt");
|
||||
swzFree(out_filename);
|
||||
printf("Constructing percentiles for forecasts - %d draws of shocks/regimes per posterior value\n",draws);
|
||||
forecast_percentile(f_out,percentiles,draws,posterior_file,p->nobs,horizon,model);
|
||||
fclose(f_out);
|
||||
}
|
||||
|
||||
if (posterior_file) fclose(posterior_file);
|
||||
|
|
|
@ -15,7 +15,7 @@ static void AddObservationFixed(PRECISION x, int *low, int *h, int *high, PRECIS
|
|||
static PRECISION Cumulative(PRECISION level, int low, int *h, PRECISION min, PRECISION max, int intervals, int sample_size);
|
||||
static PRECISION Percentile(PRECISION percentile, int low, int *h, PRECISION min, PRECISION max, int intervals, int sample_size);
|
||||
static TMatrix MakeHistogram(int low, int *h, PRECISION min, PRECISION max,int intervals, int sample_size,
|
||||
PRECISION min_out, PRECISION max_out, int bins);
|
||||
PRECISION min_out, PRECISION max_out, int bins);
|
||||
static TMatrix MakeHistogramAuto(int low, int *h, int high, PRECISION min, PRECISION max, int intervals, int sample_size, int bins);
|
||||
|
||||
|
||||
|
@ -104,14 +104,14 @@ void AddMatrixObservation(TMatrix X, TMatrixHistogram *h)
|
|||
{
|
||||
for (i=h->rows-1; i >= 0; i--)
|
||||
for (j=h->cols-1; j >= 0; j--)
|
||||
{
|
||||
h->low[i][j]=h->high[i][j]=0;
|
||||
for (k=h->intervals-1; k >= 0; k--) h->freq[i][j][k]=0;
|
||||
}
|
||||
{
|
||||
h->low[i][j]=h->high[i][j]=0;
|
||||
for (k=h->intervals-1; k >= 0; k--) h->freq[i][j][k]=0;
|
||||
}
|
||||
if (h->type == HISTOGRAM_VARIABLE)
|
||||
for (i=h->rows-1; i >= 0; i--)
|
||||
for (j=h->cols-1; j >= 0; j--)
|
||||
ElementM(h->Min,i,j)=ElementM(h->Max,i,j)=ElementM(X,i,j);
|
||||
for (j=h->cols-1; j >= 0; j--)
|
||||
ElementM(h->Min,i,j)=ElementM(h->Max,i,j)=ElementM(X,i,j);
|
||||
}
|
||||
|
||||
if (h->type == HISTOGRAM_FIXED)
|
||||
|
@ -245,12 +245,12 @@ void AddVectorObservation(TVector x, TVectorHistogram *h)
|
|||
{
|
||||
for (i=h->dim-1; i >= 0; i--)
|
||||
{
|
||||
h->low[i]=h->high[i]=0;
|
||||
for (k=h->intervals-1; k >= 0; k--) h->freq[i][k]=0;
|
||||
h->low[i]=h->high[i]=0;
|
||||
for (k=h->intervals-1; k >= 0; k--) h->freq[i][k]=0;
|
||||
}
|
||||
if (h->type == HISTOGRAM_VARIABLE)
|
||||
for (i=h->dim-1; i >= 0; i--)
|
||||
ElementV(h->Min,i)=ElementV(h->Max,i)=ElementV(x,i);
|
||||
ElementV(h->Min,i)=ElementV(h->Max,i)=ElementV(x,i);
|
||||
}
|
||||
|
||||
if (h->type == HISTOGRAM_FIXED)
|
||||
|
@ -691,57 +691,57 @@ static TMatrix MakeHistogramAuto(int low, int *h, int high, PRECISION min, PRECI
|
|||
|
||||
if (lo == intervals)
|
||||
{
|
||||
min_out=max-1.0;
|
||||
max_out=max+1.0;
|
||||
min_out=max-1.0;
|
||||
max_out=max+1.0;
|
||||
}
|
||||
else
|
||||
{
|
||||
if (high > 0)
|
||||
hi=intervals;
|
||||
else
|
||||
for (hi=intervals-1; !h[hi]; hi--);
|
||||
if (high > 0)
|
||||
hi=intervals;
|
||||
else
|
||||
for (hi=intervals-1; !h[hi]; hi--);
|
||||
|
||||
if (lo >= 0)
|
||||
if (hi < intervals)
|
||||
{
|
||||
min_out=min+lo*inc;
|
||||
max_out=min+(hi+1)*inc;
|
||||
}
|
||||
else
|
||||
{
|
||||
min_out=min+lo*inc;
|
||||
if (bins == 1)
|
||||
max_out=(1+SQRT_MACHINE_EPSILON)*max;
|
||||
else
|
||||
{
|
||||
inc=(1-SQRT_MACHINE_EPSILON)*(max - min_out)/(PRECISION)(bins-1);
|
||||
max_out=max + inc;
|
||||
}
|
||||
}
|
||||
else
|
||||
if (hi < intervals)
|
||||
{
|
||||
max_out=min+(hi+1)*inc;
|
||||
if (bins == 1)
|
||||
min_out=(1-SQRT_MACHINE_EPSILON)*min;
|
||||
else
|
||||
{
|
||||
inc=(1-SQRT_MACHINE_EPSILON)*(max_out - min)/(PRECISION)(bins-1);
|
||||
min_out=min - inc;
|
||||
}
|
||||
}
|
||||
else
|
||||
if (bins <= 2)
|
||||
{
|
||||
min_out=(1-SQRT_MACHINE_EPSILON)*min;
|
||||
max_out=(1+SQRT_MACHINE_EPSILON)*max;
|
||||
}
|
||||
else
|
||||
{
|
||||
inc=(1-SQRT_MACHINE_EPSILON)*(max_out - min)/(PRECISION)(bins-2);
|
||||
min_out=min - inc;
|
||||
max_out=max +inc;
|
||||
}
|
||||
if (lo >= 0)
|
||||
if (hi < intervals)
|
||||
{
|
||||
min_out=min+lo*inc;
|
||||
max_out=min+(hi+1)*inc;
|
||||
}
|
||||
else
|
||||
{
|
||||
min_out=min+lo*inc;
|
||||
if (bins == 1)
|
||||
max_out=(1+SQRT_MACHINE_EPSILON)*max;
|
||||
else
|
||||
{
|
||||
inc=(1-SQRT_MACHINE_EPSILON)*(max - min_out)/(PRECISION)(bins-1);
|
||||
max_out=max + inc;
|
||||
}
|
||||
}
|
||||
else
|
||||
if (hi < intervals)
|
||||
{
|
||||
max_out=min+(hi+1)*inc;
|
||||
if (bins == 1)
|
||||
min_out=(1-SQRT_MACHINE_EPSILON)*min;
|
||||
else
|
||||
{
|
||||
inc=(1-SQRT_MACHINE_EPSILON)*(max_out - min)/(PRECISION)(bins-1);
|
||||
min_out=min - inc;
|
||||
}
|
||||
}
|
||||
else
|
||||
if (bins <= 2)
|
||||
{
|
||||
min_out=(1-SQRT_MACHINE_EPSILON)*min;
|
||||
max_out=(1+SQRT_MACHINE_EPSILON)*max;
|
||||
}
|
||||
else
|
||||
{
|
||||
inc=(1-SQRT_MACHINE_EPSILON)*(max_out - min)/(PRECISION)(bins-2);
|
||||
min_out=min - inc;
|
||||
max_out=max +inc;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
|
Loading…
Reference in New Issue