Commit 39adf099 authored by Houtan Bastani's avatar Houtan Bastani
Browse files

SWZ: remove out of date code

parent 06fe0c85
sbvar_draws
sbvar_init_file
sbvar_mhm_1
sbvar_mhm_2
sbvar_probabilities
sbvar_estimation
//== Flat Independent Markov States and Simple Restrictions ==//
//-----------------------------------------------------------------------------//
//-- Read by CreateMarkovStateVariable_File() only if the passed number of --//
//-- observations is less than or equal to zero. Can be omitted if the --//
//-- passed number of observations is positive. --//
//-----------------------------------------------------------------------------//
//== Number Observations ==//
200
//== Number Independent State Variables ==//
2
//-----------------------------------------------------------------------------//
//-- state_variable[1] (1 <= i <= n_state_variables) --//
//-----------------------------------------------------------------------------//
//== Number of states for state_variable[1] ==//
3
//-----------------------------------------------------------------------------//
//-- Each column contains the parameters for a Dirichlet prior on the --//
//-- corresponding column of the transition matrix. Each element must be --//
//-- positive. For each column, the relative size of the prior elements --//
//-- determine the relative size of the elements of the transition matrix --//
//-- and overall larger sizes implies a tighter prior. --//
//-----------------------------------------------------------------------------//
//== Transition matrix prior for state_variable[1]. (n_states x n_states) ==//
10 1 1
1 10 1
1 1 10
//-----------------------------------------------------------------------------//
//-- An array of n_states integers with each entry between 1 and n_states, --//
//-- inclusive. Determines the number of quasi-free Dirichlet dimensions --//
//-- each column. Since the sum of the elements in a Dirichlet distribution --//
//-- must equal to one, the actual number of free dimensions is one less. --//
//-----------------------------------------------------------------------------//
//== Free Dirichet dimensions for state_variable[1] ==//
3 3 3
//-----------------------------------------------------------------------------//
//-- The jth restriction matrix is n_states x free[j]. Each row of the --//
//-- restriction matrix has at most one non-zero entry and the sum of each --//
//-- column of the restriction matrix must be one. If (x(1),...,x(free[j])) --//
//-- is the Dirichlet random variable for column j, then the jth column of --//
//-- the transition matrix Q is the jth restriction matrix times the --//
//-- Dirichlet random random variable. --//
//-----------------------------------------------------------------------------//
//== Column restrictions for state_variable[1] ==//
1 0 0
0 1 0
0 0 1
1 0 0
0 1 0
0 0 1
1 0 0
0 1 0
0 0 1
//-----------------------------------------------------------------------------//
//-- Allows for lagged values of the state variable to be encoded --//
//-----------------------------------------------------------------------------//
//== Number of lags encoded for state_variable[1] ==//
2
//-----------------------------------------------------------------------------//
//-- state_variable[2] --//
//-----------------------------------------------------------------------------//
//== Number of states for state_variable[2] ==//
2
//-----------------------------------------------------------------------------//
//-- Each column contains the parameters for a Dirichlet prior on the --//
//-- corresponding column of the transition matrix. Each element must be --//
//-- positive. For each column, the relative size of the prior elements --//
//-- determine the relative size of the elements of the transition matrix --//
//-- and overall larger sizes implies a tighter prior. --//
//-----------------------------------------------------------------------------//
//== Transition matrix prior for state_variable[2]. (n_states x n_states) ==//
5 1
1 5
//-----------------------------------------------------------------------------//
//-- An array of n_states integers with each entry between 1 and n_states, --//
//-- inclusive. Determines the number of quasi-free Dirichlet dimensions --//
//-- each column. Since the sum of the elements in a Dirichlet distribution --//
//-- must equal to one, the actual number of free dimensions is one less. --//
//-----------------------------------------------------------------------------//
//== Free Dirichet dimensions for state_variable[2] ==//
2 2
//-----------------------------------------------------------------------------//
//-- The jth restriction matrix is n_states x free[j]. Each row of the --//
//-- restriction matrix has at most one non-zero entry and the sum of each --//
//-- column of the restriction matrix must be one. If (x(1),...,x(free[j])) --//
//-- is the Dirichlet random variable for column j, then the jth column of --//
//-- the transition matrix Q is the jth restriction matrix times the --//
//-- Dirichlet random random variable. --//
//-----------------------------------------------------------------------------//
//== Column restrictions for state_variable[2] ==//
1 0
0 1
1 0
0 1
//-----------------------------------------------------------------------------//
//-- Allows for lagged values of the state variable to be encoded --//
//-----------------------------------------------------------------------------//
//== Number of lags encoded for state_variable[2] ==//
0
This diff is collapsed.
This diff is collapsed.
#include "switch_opt.h"
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "modify_for_mex.h"
/* //====== Static Global Variables ====== ansi-c*/
static TStateModel *Model=(TStateModel*)NULL;
static PRECISION *buffer=(PRECISION*)NULL;
static PRECISION *ModifiedFreeParameters=(PRECISION*)NULL;
static PRECISION *FreeParameters_Q=(PRECISION*)NULL;
static int NumberFreeParameters_Q=0;
static PRECISION *FreeParameters_Theta=(PRECISION*)NULL;
static int NumberFreeParameters_Theta=0;
void SetupObjectiveFunction(TStateModel *model, PRECISION *Modified, PRECISION *FreeQ, PRECISION *FreeTheta)
{
if (buffer) swzFree(buffer);
Model=model;
FreeParameters_Q=FreeQ;
NumberFreeParameters_Q=NumberFreeParametersQ(model);
FreeParameters_Theta=FreeTheta;
NumberFreeParameters_Theta=model->routines->pNumberFreeParametersTheta(model);
ModifiedFreeParameters=Modified;
}
void SetupObjectiveFunction_new(TStateModel *model, int FreeTheta_Idx, int FreeQ_Idx, int Modified_Idx)
{
if (buffer) swzFree(buffer);
Model=model;
NumberFreeParameters_Q=NumberFreeParametersQ(model);
NumberFreeParameters_Theta=model->routines->pNumberFreeParametersTheta(model);
buffer=(PRECISION*)swzMalloc((NumberFreeParameters_Q + NumberFreeParameters_Theta)*sizeof(PRECISION));
FreeParameters_Q=buffer+FreeQ_Idx;
FreeParameters_Theta=buffer+FreeTheta_Idx;
ModifiedFreeParameters=buffer+Modified_Idx;
}
PRECISION PosteriorObjectiveFunction(PRECISION *x, int n)
{
if (x != ModifiedFreeParameters) memmove(ModifiedFreeParameters,x,n*sizeof(PRECISION));
ConvertFreeParametersToQ(Model,FreeParameters_Q);
ConvertFreeParametersToTheta(Model,FreeParameters_Theta);
return -LogPosterior_StatesIntegratedOut(Model);
/* //PRECISION lp_Q, lp_Theta, li; ansi-c*/
/* //FILE *f_out; ansi-c*/
/* //lp_Q=LogPrior_Q(Model); ansi-c*/
/* //lp_Theta=LogPrior_Theta(Model); ansi-c*/
/* //li=LogLikelihood_StatesIntegratedOut(Model); ansi-c*/
/* //if (isnan(lp_Q) || isnan(lp_Theta) || isnan(li)) ansi-c*/
/* // { ansi-c*/
/* // f_out=fopen("tmp.tmp","wt"); ansi-c*/
/* // Write_VAR_Specification(f_out,(char*)NULL,Model); ansi-c*/
/* // WriteTransitionMatrices(f_out,(char*)NULL,"Error: ",Model); ansi-c*/
/* // Write_VAR_Parameters(f_out,(char*)NULL,"Error: ",Model); ansi-c*/
/* // fprintf(f_out,"LogPrior_Theta(): %le\n",lp_Theta); ansi-c*/
/* // fprintf(f_out,"LogPrior_Q(): %le\n",lp_Q); ansi-c*/
/* // fprintf(f_out,"LogLikelihood_StatesIntegratedOut(): %le\n",li); ansi-c*/
/* // fprintf(f_out,"Posterior: %le\n\n",lp_Q+lp_Theta+li); ansi-c*/
/* // fclose(f_out); ansi-c*/
/* // swzExit(0); ansi-c*/
/* // } ansi-c*/
/* //return -(lp_Q+lp_Theta+li); ansi-c*/
}
PRECISION PosteriorObjectiveFunction_csminwel(double *x, int n, double **args, int *dims)
{
return PosteriorObjectiveFunction(x,n);
}
void PosteriorObjectiveFunction_npsol(int *mode, int *n, double *x, double *f, double *g, int *nstate)
{
*f=PosteriorObjectiveFunction(x,*n);
}
PRECISION MLEObjectiveFunction(PRECISION *x, int n)
{
if (x != ModifiedFreeParameters) memmove(ModifiedFreeParameters,x,n*sizeof(PRECISION));
ConvertFreeParametersToQ(Model,FreeParameters_Q);
ConvertFreeParametersToTheta(Model,FreeParameters_Theta);
return -LogLikelihood_StatesIntegratedOut(Model);
}
PRECISION MLEObjectiveFunction_csminwel(double *x, int n, double **args, int *dims)
{
return MLEObjectiveFunction(x,n);
}
void MLEObjectiveFunction_npsol(int *mode, int *n, double *x, double *f, double *g, int *nstate)
{
*f=MLEObjectiveFunction(x,*n);
}
PRECISION MLEObjectiveFunction_LogQ(PRECISION *x, int n)
{
if (x != ModifiedFreeParameters) memmove(ModifiedFreeParameters,x,n*sizeof(PRECISION));
ConvertLogFreeParametersToQ(Model,FreeParameters_Q);
ConvertFreeParametersToTheta(Model,FreeParameters_Theta);
return -LogLikelihood_StatesIntegratedOut(Model);
}
PRECISION MLEObjectiveFunction_LogQ_csminwel(double *x, int n, double **args, int *dims)
{
return MLEObjectiveFunction_LogQ(x,n);
}
#ifndef __MARKOV_SWITCHING_OPTIMIZATION__
#define __MARKOV_SWITCHING_OPTIMIZATION__
#include "switch.h"
void SetupObjectiveFunction(TStateModel *model, PRECISION *MFPparms, PRECISION *FTMparms, PRECISION *FMSparms);
PRECISION PosteriorObjectiveFunction(PRECISION *x, int n);
PRECISION PosteriorObjectiveFunction_csminwel(double *x, int n, double **args, int *dims);
void PosteriorObjectiveFunction_npsol(int *mode, int *n, double *x, double *f, double *g, int *nstate);
PRECISION MLEObjectiveFunction(PRECISION *x, int n);
PRECISION MLEObjectiveFunction_csminwel(double *x, int n, double **args, int *dims);
void MLEObjectiveFunction_npsol(int *mode, int *n, double *x, double *f, double *g, int *nstate);
PRECISION MLEObjectiveFunction_LogQ(PRECISION *x, int n);
PRECISION MLEObjectiveFunction_LogQ_csminwel(double *x, int n, double **args, int *dims);
#endif
This diff is collapsed.
#include "switch.h"
/*
Base routines for reading/writing Markov state variables and transition
matrices in native ascii format.
*/
TMarkovStateVariable* ReadMarkovSpecification_SV(FILE *f_in, char *idstring, int nobs);
int WriteMarkovSpecification_SV(FILE *f_out, TMarkovStateVariable *sv, char *idstring);
int ReadTransitionMatrices_SV(FILE *f_in, TMarkovStateVariable* sv, char *header, char *idstring);
int WriteTransitionMatrices_SV(FILE *f_out, TMarkovStateVariable* sv, char *header, char *idstring);
int ReadBaseTransitionMatrices_SV(FILE *f_out, TMarkovStateVariable *sv, char *header, char *idstring);
int WriteBaseTransitionMatrices_SV(FILE *f_out, TMarkovStateVariable *sv, char *header, char *idstring);
int ReadBaseTransitionMatricesFlat_SV(FILE *f_out, TMarkovStateVariable *sv);
int WriteBaseTransitionMatricesFlat_SV(FILE *f_out, TMarkovStateVariable *sv, char *fmt);
void WriteBaseTransitionMatricesFlat_Headers_SV(FILE *f_out, TMarkovStateVariable* sv, char *idstring);
/*
Routines for reading/writing Markov state variables and transition matrices
from TStateModel. Calls base routines.
*/
TMarkovStateVariable* ReadMarkovSpecification(FILE *f, char *filename);
int WriteMarkovSpecification(FILE *f, char *filename, TStateModel *model);
int ReadTransitionMatrices(FILE *f, char *filename, char *header, TStateModel *model);
int WriteTransitionMatrices(FILE *f, char *filename, char *header, TStateModel *model);
int ReadStates(FILE *f, char *filename, char *header, TStateModel *model);
int WriteStates(FILE *f, char *filename, char *header, TStateModel *model);
int ReadBaseTransitionMatrices(FILE *f, char *filename, char *header, TStateModel *model);
int WriteBaseTransitionMatrices(FILE *f, char *filename, char *header, TStateModel *model);
int ReadBaseTransitionMatricesFlat(FILE *f, TStateModel *model);
int WriteBaseTransitionMatricesFlat(FILE *f, TStateModel *model, char *fmt);
/*
Read flat markov state variable specification from file.
*/
TMarkovStateVariable* CreateMarkovStateVariable_File(FILE *f, char *filename, int nobs);
#include "swzmatrix.h"
#include "dw_rand.h"
#include "dw_parse_cmd.h"
#include "dw_ascii.h"
#include "VARbase.h"
#include "VARio.h"
#include "switch.h"
#include "switchio.h"
#include "command_line_VAR.h"
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <time.h>
#include "modify_for_mex.h"
int main(int nargs, char **args)
{
TStateModel *model;
T_VAR_Parameters *p;
FILE *f_out;
char *filename;
int count, begin_time, end_time, tuning, burn_in, iterations, check, period=1000, seed, output, thinning,
nd1;
TVARCommandLine *cmd=(TVARCommandLine*)NULL;
char *include_help[]={"-di","-do","-fs","-fp","-ph","-MLE",(char*)NULL},
*additional_help[]={
"-ft <tag>",
"Tag for input file. Input file name is est_final_<tag>.dat.",
"-fto <tag>",
"Tag for output file. Output file names are draws_<tag>.dat and headers_<tag>.dat. Default is -ft <tag>.",
"-mh <integer>",
"Tuning period for Metropolis-Hasting draws (default value = 30000)",
"-b <integer>",
"Burn-in period (default value = 0.1 * (number of iterations))",
"-i <integer>",
"Number of draw (default value = 1000), number saved is (-i)/(-t)",
"-t <integer>",
"Thinning factor. Only one in t draws are written to file.",
"-nd1",
"Normalize diagonal of A0 to one (flat output only)",
"-gs <integer>",
"Seed value for generator - 0 gets seed from clock (default value = 0)",
(char*)NULL,
(char*)NULL};
/* //=== Help Screen === ansi-c*/
if (dw_FindArgument_String(nargs,args,"h") != -1)
{
printf("print_draws <options>\n");
PrintHelpMessages(stdout,include_help,additional_help);
return 0;
}
/* //=== Get seed, tuning peroid, burn-in period, number of iterations, and thinning factor ansi-c*/
constant_seed=dw_ParseInteger_String(nargs,args,"cseed",0);
seed=dw_ParseInteger_String(nargs,args,"gs",0);
tuning=dw_ParseInteger_String(nargs,args,"mh",30000);
iterations=dw_ParseInteger_String(nargs,args,"i",1000);
burn_in=dw_ParseInteger_String(nargs,args,"b",iterations/10);
thinning=dw_ParseInteger_String(nargs,args,"t",1);
nd1=(dw_FindArgument_String(nargs,args,"nd1") >= 0) ? 1 : 0;
/* //=== Initialize random number generator ansi-c*/
dw_initialize_generator(seed);
/* //=== Setup model and initial parameters ansi-c*/
printf("Reading data...\n");
if (!(model=CreateTStateModelFromEstimateFinal(nargs,args,&cmd)))
{
swz_fprintf_err("Unable to read model or parameters\n");
swzExit(1);
}
p=(T_VAR_Parameters*)(model->theta);
/* //=== Open header file and print headers ansi-c*/
filename=CreateFilenameFromTag("%sheader_%s.dat",cmd->out_tag,cmd->out_directory);
f_out=fopen(filename,"wt");
swzFree(filename);
WriteBaseTransitionMatricesFlat_Headers_SV(f_out,model->sv,"");
Write_VAR_ParametersFlat_Headers(f_out,model);
fprintf(f_out,"\n");
fclose(f_out);
/* //=== Open output file ansi-c*/
filename=CreateFilenameFromTag("%sdraws_%s.dat",cmd->out_tag,cmd->out_directory);
f_out=fopen(filename,"wt");
swzFree(filename);
/* // Burn-in period with calibration of jumping parameters ansi-c*/
printf("Calibrating jumping parameters - %d draws\n",tuning);
begin_time=(int)time((time_t*)NULL);
AdaptiveMetropolisScale(model,tuning,1000,1,(FILE*)NULL); /* tuning iterations - 1000 iterations before updating - verbose ansi-c*/
end_time=(int)time((time_t*)NULL);
printf("Elapsed Time: %d seconds\n",end_time - begin_time);
/* // Reset parametrers ansi-c*/
if (!ReadTransitionMatrices((FILE*)NULL,cmd->parameters_filename_actual,cmd->parameters_header_actual,model)
|| !Read_VAR_Parameters((FILE*)NULL,cmd->parameters_filename_actual,cmd->parameters_header_actual,model))
printf("Unable to reset parameters after tuning\n");
/* // Burn-in period ansi-c*/
printf("Burn-in period - %d draws\n",burn_in);
for (check=period, count=1; count <= burn_in; count++)
{
DrawAll(model);
if (count == check)
{
check+=period;
printf("%d iterations completed out of %d\n",count,burn_in);
}
}
end_time=(int)time((time_t*)NULL);
printf("Elapsed Time: %d seconds\n",end_time - begin_time);
ResetMetropolisInformation(p);
/* // Simulation ansi-c*/
printf("Simulating - %d draws\n",iterations);
for (check=period, output=thinning, count=1; count <= iterations; count++)
{
DrawAll(model);
if (count == output)
{
WriteBaseTransitionMatricesFlat(f_out,model,"%lf ");
if (nd1)
Write_VAR_ParametersFlat_A0_Diagonal_One(f_out,model,"%lf ");
else
Write_VAR_ParametersFlat(f_out,model,"%lf ");
fprintf(f_out,"\n");
output+=thinning;
}
if (count == check)
{
check+=period;
printf("%d(%d) iterations completed out of %d(%d)\n",count,thinning,iterations,thinning);
}
}
end_time=(int)time((time_t*)NULL);
printf("Elapsed Time: %d seconds\n",end_time - begin_time);
/* // clean up ansi-c*/
fclose(f_out);
FreeStateModel(model);
Free_VARCommandLine(cmd);
}
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
#ifndef __VAR_INPUT_OUTPUT__
#define __VAR_INPUT_OUTPUT__
#include "switch.h"
#include "VARbase.h"
void Write_VAR_Specification(FILE *f, char *filename, TStateModel *model);
TStateModel* Read_VAR_Specification(FILE *f, char *filename);
int Write_VAR_Parameters(FILE *f, char *filename, char *id, TStateModel *model);
int Read_VAR_Parameters(FILE *f, char *filename, char *id, TStateModel *model);
int Read_VAR_ParametersFlat(FILE *f_in, TStateModel *model);
int Write_VAR_ParametersFlat(FILE *f, TStateModel *model, char *fmt);
int Write_VAR_ParametersFlat_Headers(FILE *f_out, TStateModel *model);
int Write_VAR_ParametersFlat_A0_Diagonal_One(FILE *f, TStateModel *model, char *fmt);
void ReadAllParameters(FILE *f, char *filename, char *id, TStateModel *model);
void WriteAllParameters(FILE *f, char *filename, char *id, TStateModel *model);
/* //T_VAR_Parameters* Create_VAR_Parameters_File(FILE *f, char *filename, TMarkovStateVariable *sv); ansi-c*/
/* //TStateModel* CreateStateModel_VAR_File(FILE *f, char *filename); ansi-c*/
/* //void PrintParametersVAR(FILE *f_out, TStateModel *model); ansi-c*/
/* //void Write_VAR_Info(FILE *f, char *filename, T_VAR_Parameters *p); ansi-c*/
#endif
This diff is collapsed.
#ifndef __VAR_INPUT_OUTPUT_MATLAB__
#define __VAR_INPUT_OUTPUT_MATLAB__
#include "switch.h"
#include "VARbase.h"
TStateModel* Combine_matlab_standard(char *inmatlab, char *instandard);
void ReadConstantParameters(char *filename, TStateModel *model);
TStateModel* CreateStateModel_VAR_matlab(char *filename);
#endif
This diff is collapsed.
#ifndef __command_line_VAR
#define __command_line_VAR
#include "switch.h"
#include <stdio.h>
char* CreateFilenameFromTag(char *fmt, char *tag, char *dir);
char* CreatePath(char *path);
void PrintHelpMessages(FILE *f, char **include, char **additional);
typedef struct
{
char *in_directory; /* -di ansi-c*/
char *in_tag; /* -ft ansi-c*/
char *specification_filename; /* -fs ansi-c*/
char *parameters_filename; /* -fp ansi-c*/
char *parameters_header; /* -ph ansi-c*/
char *specification_filename_actual;
char *parameters_filename_actual;
char *parameters_header_actual;
char *out_directory; /* -do ansi-c*/
char *out_tag; /* -fto (default from -ft) ansi-c*/
char *out_header; /* -pho (default from -ph) ansi-c*/
} TVARCommandLine;
TVARCommandLine* Create_VARCommandLine(void);
void Free_VARCommandLine(TVARCommandLine *cmd);
TVARCommandLine* Base_VARCommandLine(int nargs, char **args, TVARCommandLine *cmd);
void EstimateFinal_VARCommandLine_Help(FILE *f);
TStateModel* CreateTStateModelFromEstimateFinal(int nargs, char **args, TVARCommandLine **p_cmd);
TStateModel* CreateTStateModelForEstimate(int nargs, char **args, TVARCommandLine **p_cmd);
#endif
#include "switch.h"
#include "switchio.h"
#include "VARbase.h"
#include "VARio.h"
#include "VARio_matlab.h"
#include "dw_error.h"
#include "dw_ascii.h"
#include <stdlib.h>
#include <string.h>
#include "modify_for_mex.h"
/*
Creates a standard initialization file from the matlab and specification file.
*/
int main(int nargs, char **args)
{
TStateModel *model;
FILE *f_out;
char *filename, *fmt="init_%s.dat", *header="Initial: ";
dw_SetTerminalErrors(ALL_ERRORS);
dw_SetVerboseErrors(ALL_ERRORS);
if (nargs != 4)
{
swz_fprintf_err("Syntax:\n create_init_file <matlab filename> <specs filename> <file tag>\n");
swzExit(0);
}
model=Combine_matlab_standard(args[1],args[2]);
ReadConstantParameters(args[1],model);