Thank you for the advice. is now correct?TomDoan wrote:Multiplying the log by the inflation factor is *really* not the right thing to do. The first part of that is actually computing v^(-.5) (using exp(-.5*log(v)) so you either divide the whole factor by sqrt(1+alfa*dummy) or add log(1+alfa*dummy) to the v before it gets multiplied by the -.5.
Code: Select all
compute %MSVARFVec(i)=exp(-.5*MSVARLogDetV(state)/sqrt(1+alfa*dummy(t))-.5*%qform(MSVARSigmaInv(state)/(1+alfa*dummy(t)),u)
Thank you for pointing out the mistake.TomDoan wrote:No. You're dividing in the first term inside the exp, not outside.
Code: Select all
compute %MSVARFVec(i)=exp(-.5*(MSVARLogDetV(state)+Log(1+alfa*dummy(t)))-.5*%qform(MSVARSigmaInv(state)/(1+alfa*dummy(t)),u)
Code: Select all
cal(m) 1948
open data filardo.dat
data(format=prn,org=cols) 1948:1 1992:8 year month ip cli dfi xli sp fed tspread
*
set ipgr = 100.0*log(ip/ip{1})
*
set g = ipgr
*
* Transform the information variables
*
set cligr = 100.0*log(cli/cli{1})
set spgr = 100.0*log(sp/sp{1})
set xlidiff = xli-xli{1}
set ffdiff = fed-fed{1}
set tspdiff = tspread-tspread{1}
set dfismooth = dfi+2*dfi{1}+2*dfi{2}+dfi{3}
*
* And center them to zero means
*
dofor s = cligr spgr xlidiff ffdiff tspdiff dfismooth
diff(center) s
end dofor s
*******************************************************
*
* Common setup
*
compute nlags=4
*
source msvarsetup.src
compute gstart=1948:6,gend=1992:8
@MSVARSetup(lags=nlags)
# g
*********************************************************************
*
*
nonlin(parmset=common) mu phi sigma alfa
@MSVARinitial gstart gend
*********************************************************************
*
* Time varying probabilities
*
dec equation p1eq p2eq
dec vector v1 v2
nonlin(parmset=tv) v1 v2
*********************************************************************
*
* This overrides the standard fixed definition
*
function %MSVARPmat time
type rect %MSVARPmat
type int time
*
local integer i j
local rect pexpand
local real z
*
compute p(1,1)=%(z=exp(%dot(%eqnxvector(p1eq,time),v1)),z/(1+z))
compute p(1,2)=%(z=exp(%dot(%eqnxvector(p2eq,time),v2)),1/(1+z))
compute %MSVARInitTransition()
*
if nexpand==nstates {
dim pexpand(nstates,nstates-1)
ewise pexpand(i,j)=%if(i==nstates,-1.0,i==j)
compute %MSVARPmat=pexpand*p
ewise %MSVARPmat(i,j)=%if(i==nstates,1.0+%MSVARPmat(i,j),%MSVARPmat(i,j))
}
else {
dim %MSVARPmat(nexpand,nexpand)
ewise %MSVARPmat(i,j)=MSVARTransProbs(MSVARTransLookup(i,j))
}
end
*********************************************************************
*
* Initialize the probabilities using only the intercepts
*
function TVPInit
type vector TVPInit
*
local rect a
local integer i j
*
compute p=||%logistic(v1(1),1.0),1-%logistic(v2(1),1.0)||
compute TVPInit=%MSVARInit()
end
*********************************************************************
set dummy = t<=1959:12
dec real alfa
*
* Override the %MSVARFVec function. This returns the expanded set of
* likelihood functions across all state combinations.
*
function %MSVARFVec time
type vector %MSVARFVec
type integer time
*
local integer i j state
local vect u
local real iconst
*
dim %MSVARFVec(nexpand)
*
compute iconst=(2*%pi)^(-.5*nvar)
do i=1,nexpand
compute state=MSVARLagState(i,1)
compute u=%xt(y,time)-mu(state)(1)
do j=1,nlags
compute statej=MSVARLagState(i,j+1),$
u=(u-phi(j)*(%xt(y,time-j)-mu(statej)(1))
end do j
compute state=MSVARLagState(i,1)
compute %MSVARFVec(i)=iconst*exp(-.5*(MSVARLogDetV(state)+Log(1+alfa*dummy(t)))-.5*%qform(MSVARSigmaInv(state)/(1+alfa*dummy(t)),u)
end do i
*
end
frml msvarf = %MSVARPMat(t),log(%MSVarProb(t))
*********************************************************************
*
* Define the logistic index for the transitions
*
equation p1eq *
# constant cligr{1}
equation p2eq *
# constant cligr{1}
*
* Standard guess values aiming to make state 1 the low growth state.
*
compute v1=log(.7/.3)~~0.0
compute v2=log(.95/.05)~~0.0
compute mu(1)=-1.7,mu(2)=.3
*
maximize(parmset=tv+common,start=(pstar=TVPInit()),$
method=bhhh,iters=20) logltvtp gstart gend
Thank you so much for your help.TomDoan wrote:That looks fine. You can do an obvious check by pegging alfa to zero and seeing whether you get the same results as you do for the model without the heteroscedasticity correction.
Dear Tom, sorry for posting question again. After applying the code to my own exercise, I got an error message says: ## MAT13. Store into Out-of-Range Matrix or Series ElementTomDoan wrote:That looks fine. You can do an obvious check by pegging alfa to zero and seeing whether you get the same results as you do for the model without the heteroscedasticity correction.
Code: Select all
CALENDAR(Q) 1875:1
DATA(FORMAT=XLSX,ORG=COLUMNS) 1875:01 2010:04 quarter recession rgdp rep yr34
set demo = 1-rep
set yr12 = 1-yr34
set grgdp = 100.0*log(rgdp/rgdp{1})
graph(footer="Figure 1. Log Growth Rate of Quarter Real GDP")
# grgdp
states(smpl=yr34) grgdp
@nbercycles(up=ups,down=downs)
states(smpl=rep) ups 1950:01 *
set g = grgdp
sstats(mean,smpl=ups{1}) 1950:01 * ups>>pups
display pups
sstats(mean,smpl=downs{1}) 1950:01 * downs>>pdowns
display pdowns
stats(smpl=downs) g /
*******************************************************
*
* Common setup
*
compute nlags=4
*
source msvarsetup.src
compute gstart=1876:02,gend=2010:04
@MSVARSetup(lags=nlags,switch=m)
# g
*********************************************************************
nonlin(parmset=common) mu phi sigma alfa
@MSVARinitial gstart gend
************************************************************************
*
* Time varying probabilities
*
dec equation p1eq p2eq
dec vector v1 v2
nonlin(parmset=tv) v1 v2
*********************************************************************
*
* This overrides the standard fixed definition
*
function %MSVARPmat time
type rect %MSVARPmat
type int time
*
local integer i j
local rect pexpand
local real z
*
compute p(1,1)=%(z=exp(%dot(%eqnxvector(p1eq,time),v1)),z/(1+z))
compute p(1,2)=%(z=exp(%dot(%eqnxvector(p2eq,time),v2)),1/(1+z))
compute %MSVARInitTransition()
*
if nexpand==nstates {
dim pexpand(nstates,nstates-1)
ewise pexpand(i,j)=%if(i==nstates,-1.0,i==j)
compute %MSVARPmat=pexpand*p
ewise %MSVARPmat(i,j)=%if(i==nstates,1.0+%MSVARPmat(i,j),%MSVARPmat(i,j))
}
else {
dim %MSVARPmat(nexpand,nexpand)
ewise %MSVARPmat(i,j)=MSVARTransProbs(MSVARTransLookup(i,j))
}
end
*********************************************************************
*
* Initialize the probabilities using only the intercepts
*
function TVPInit
type vector TVPInit
*
local rect a
local integer i j
*
compute p=||%logistic(v1(1),1.0),1-%logistic(v2(1),1.0)||
compute TVPInit=%MSVARInit()
end
*********************************************************************
set dummy = t<=1949:04
dec real alfa
*
* Override the %MSVARFVec function. This returns the expanded set of
* likelihood functions across all state combinations.
*
function %MSVARFVec time
type vector %MSVARFVec
type integer time
*
local integer i j state
local vect u
local real iconst
*
compute iconst=(2*%pi)^(-.5*nvar)
do i=1,nexpand
compute state=MSVARLagState(i,1)
compute u=%xt(y,time)-mu(state)(1)
do j=1,nlags
compute statej=MSVARLagState(i,j+1),$
u=u-phi(j)*(%xt(y,time-j)-mu(statej)(1))
end do j
compute state=MSVARLagState(i,1)
compute %MSVARFVec(i)=iconst*exp(-.5*(MSVARLogDetV(state)+Log(1+alfa*dummy(t)))-.5*%qform(MSVARSigmaInv(state)/(1+alfa*dummy(t)),u))
end do i
*
end
frml msvarf = %MSVARPMat(t),log(%MSVarProb(t))
*********************************************************************
*
* Define the logistic index for the transitions
*
equation p1eq *
# yr34 yr12
equation p2eq *
# yr34 yr12
*
sstats(mean,smpl=yr12.and.downs{1}) / downs>>pdowns
display pdowns
sstats(mean,smpl=yr12.and.ups{1}) / ups>>pups
display pups
* Standard guess values aiming to make state 1 the low growth state.
*
compute v1=log(.8/.2)~~log(.85/.15)
compute v2=log(.94/.06)~~log(.9/.1)
compute mu(1)=-0.8,mu(2)=1.5
*
maximize(parmset=tv+common,start=(pstar=TVPInit()),$
method=bhhh,iters=400,pmethod=simplex,piters=20) msvarf gstart gend
* Do graph of NBER recessions vs the probability of slow growth state
*
@NBERCycles(down=recession)
@MSVARsmoothed gstart gend psmooth
set pdown gstart gend = psmooth(t)(1)
graph(shade=recession)
# pdown
Thank you so much for pointing out the missing instruction for me. The code now works but the estimates on transition indices (V1(1), V2(1)) are surprisingly large. I guess this may be due to the initial values used for the estimation. I tried several initial values but none can fix the problem. Could you please kind take a look at my data (attached) and let me know how you will pick up the initial values? My data cover the period from 1875:Q1 to 2010:Q4, and based on the recessions and expansions defined by NBER, the average growth rate for recession period is -0.82 and is 2.1 for expansion,before WW2, up to 1949:Q4. while, the average growth rate for recession period is -0.54 and is 1.0 for expansion for the after-war period. It seems the gap between expansion and recession is getting smaller in the later period. Is it necessary to account this change into the estimation? Thank you in advanceTomDoan wrote:You never dimension %MSVARFVEC inside the function.
Code: Select all
OPEN DATA "C:\Users\Sam\Desktop\historyofrealgdp.xlsx"
CALENDAR(Q) 1875:1
DATA(FORMAT=XLSX,ORG=COLUMNS) 1875:01 2010:04 quarter recession rgdp rep yr34
set demo = 1-rep
set yr12 = 1-yr34
set grgdp = 100.0*log(rgdp/rgdp{1})
graph(footer="Figure 1. Log Growth Rate of Quarter Real GDP")
# grgdp
states(smpl=recession) grgdp 1950:01 *
states(smpl=.not.recession) grgdp 1950:01 *
@nbercycles(up=ups,down=downs)
states(smpl=rep) ups 1950:01 *
set g = grgdp
sstats(mean,smpl=ups{1}) 1950:01 * ups>>pups
display pups
sstats(mean,smpl=downs{1}) 1950:01 * downs>>pdowns
display pdowns
stats(smpl=downs) g /
*******************************************************
*
* Common setup
*
compute nlags=4
*
source msvarsetup.src
compute gstart=1876:02,gend=2010:04
@MSVARSetup(lags=nlags,switch=m)
# g
*********************************************************************
nonlin(parmset=common) mu phi sigma alfa mustar
@MSVARinitial gstart gend
************************************************************************
*
* Time varying probabilities
*
dec equation p1eq p2eq
dec vector v1 v2
nonlin(parmset=tv) v1 v2
*********************************************************************
*
* This overrides the standard fixed definition
*
function %MSVARPmat time
type rect %MSVARPmat
type int time
*
local integer i j
local rect pexpand
local real z
*
compute p(1,1)=%(z=exp(%dot(%eqnxvector(p1eq,time),v1)),z/(1+z))
compute p(1,2)=%(z=exp(%dot(%eqnxvector(p2eq,time),v2)),1/(1+z))
compute %MSVARInitTransition()
*
if nexpand==nstates {
dim pexpand(nstates,nstates-1)
ewise pexpand(i,j)=%if(i==nstates,-1.0,i==j)
compute %MSVARPmat=pexpand*p
ewise %MSVARPmat(i,j)=%if(i==nstates,1.0+%MSVARPmat(i,j),%MSVARPmat(i,j))
}
else {
dim %MSVARPmat(nexpand,nexpand)
ewise %MSVARPmat(i,j)=MSVARTransProbs(MSVARTransLookup(i,j))
}
end
*********************************************************************
*
* Initialize the probabilities using only the intercepts
*
function TVPInit
type vector TVPInit
*
local rect a
local integer i j
*
compute p=||%logistic(v1(1),1.0),1-%logistic(v2(1),1.0)||
compute TVPInit=%MSVARInit()
end
*********************************************************************
set dummy = t<=1949:04
dec real alfa
dec vect mustar(2)
*
* Override the %MSVARFVec function. This returns the expanded set of
* likelihood functions across all state combinations.
*
function %MSVARFVec time
type vector %MSVARFVec
type integer time
*
local integer i j state
local vect u
local real iconst
*
dim %MSVARFVec(nexpand)
*
compute iconst=(2*%pi)^(-.5*nvar)
do i=1,nexpand
compute state=MSVARLagState(i,1)
compute u=%xt(y,time)-mu(state)(1)-mustar(state)*dummy(t)
do j=1,nlags
compute statej=MSVARLagState(i,j+1),$
u=u-phi(j)*(%xt(y,time-j)-mu(statej)(1)-mustar(statej)*dummy(time-j))
end do j
compute state=MSVARLagState(i,1)
compute %MSVARFVec(i)=iconst*exp(-.5*(MSVARLogDetV(state)+Log(1+alfa*dummy(t)))-.5*%qform(MSVARSigmaInv(state)/(1+alfa*dummy(t)),u))
end do i
*
end
********************************************************************************
frml msvarf = %MSVARPMat(t),log(%MSVarProb(t))
*********************************************************************
*
* Define the logistic index for the transitions
*
equation p1eq *
# yr34 yr12
equation p2eq *
# yr34 yr12
*
sstats(mean,smpl=yr12.and.downs{1}) / downs>>pdowns
display pdowns
sstats(mean,smpl=yr12.and.ups{1}) / ups>>pups
display pups
* Standard guess values aiming to make state 1 the low growth state.
*
compute v1=log(.8/.2)~~log(.85/.15)
compute v2=log(.94/.06)~~log(.9/.1)
compute mu(1)=-0.8, mustar(1)=0.0
compute mu(2)=1.5, mustar(2)=0.0
*
maximize(parmset=tv+common,start=(pstar=TVPInit()),$
method=bhhh,iters=100,pmethod=simplex,piters=5) msvarf gstart gend
* Do graph of NBER recessions vs the probability of slow growth state
*
@NBERCycles(down=recession)
@MSVARsmoothed gstart gend psmooth
set pdown gstart gend = psmooth(t)(1)
graph(shade=recession)
# pdown
Thank you for your quick reply. I forgot to delete the MUSTAR I used in a different code to catch the narrowing gap between expansion and recession after WW2. Following your suggestion, I changed (T) to (TIME) in %MSVARFVEC function and gave a guess value (3.0) to ALFA. Unfortunately, the results are still strange. I have quarterly data from more than 100 years (1857:Q1 to 2010:Q4), and there are much more ups and downs in the early part of the sample period. I think that makes more difficult to describe the data through how to pick up the guess values and number of lags in autocorrelation part as the economy may behave very differently in first 50 years and remaining years. Could you please kindly let me know how you will solve the problems if you were having such data set? Many ThanksTomDoan wrote:There are all kinds of problems with what you posted. First, you have MUSTAR added to the parameter set before you've even declared it as a VECTOR. ALFA isn't given a guess value. But the main problem is that you are using (T) rather than (TIME) several times in your %MSVARFVEC function.
Thank you for the detailed reply. can I ask one more question? In most cases, the variables used in the probability transition equation are exogenous. I am wondering whether endogenous variables can be used in the transition equation. In most political economic studies, authors seem to ignore the endogeneity between politics and economics. I would like to incorporate this endogeneity into my study. For example, in UK, government can call election based on recent past economic condition; and the recent past economic condition can influence the outcome of the election. In other words, the political variables used in the transition equation can be determined by the probability of the economic status in the recent past. Could you please kindly share your suggestion on how I can modify the codes to incorporate the endogeneity into my estimation? Many Thanks in advanceTomDoan wrote:First of all, you probably do want the MUSTAR's in there, as the ups and downs are clearly higher in the early part of the data set. (Simply correcting for heteroscedasticity won't change the fact that your model switches the mean and the means are different as well).
This model tends to pick up a relatively few really extreme down periods. If those all turn out to be in "year12" periods, then the transition probability to the low regime in the "year34" can be zero. To create zero (or one) in a logistic requires "infinity" index values. Even if there are some data points that are borderline regime 1 that aren't in "year12", the likelihood might be improved overall by effectively ignoring those: confusion as to the regime can hurt the likelihood, so it may be better for the likelihood to be certain and be "wrong".