Parameter stability testing for Austria
Show older comments
Dear all,
I am having a major issue at testing for parameter stability, using ESS procedure, developed by Inoue and Rossi. For Slovenia and Italy, I obtain perfect result, a positive SupLR value, that is result of Andrews (1993) QLR test. However, for Austria, the value is negative! If I make a modification to the formula for Chow test, than the SupLR value becomes positive, however, I can't change the methodology for consistency reasons.
Does anyone have any idea what is wrong with Austria?
Thank you so much!
Here is the formula, line 53 in attached file:
- original with minus in the beginning:
chow=-2* (L-(L1*(t2/T)+(1-t2/T)*L2));
2. modified without minus in the beginning:
chow=2* (L-(L1*(t2/T)+(1-t2/T)*L2));
Answers (1)
No way anybody can tell other than by going back to the beginning for the specific data set.
But, the problem is that the quantity
(L1*(t2/T)+(1-t2/T)*L2)
is less than L while it must be greater than L in order for -2 times the difference to be positive.
We, of course, have no idea of what any of those terms are nor how they're obtained but something went awry long before now it would appear.
Although the t1,12 over T terms are weighting L1, L2 by a fractional time duration over the total; possibly a not-too-bad problem might be the time frames provided aren't consistent? That should be relatively easy to check and if were to be a problem, relatively easy fix. More than that that there's something bum about the L values is probably a lot more trouble to debug.
The above would be somewhat simpler to read if the terms were written consistently as
chow=-2*(L-(L1*(t2/T) + L2*(1-t2/T)));
NOTA BENE that t2/T + (1-t2/T) --> 1 + t2/T -t2/T --> 1
Consequently, for the end result to be positive we've already noted that the whole expression in L must be negative. That also then implies that (L1 + L2) must be > L.
This is undoubtedly more iffy code full of all the previously identified issues, but let's see what it looks like.
It would have been MUCH simpler/better to have provided the code where could be read, btw...
L=readlines('test.m'); % load the m-file
L=L(strlength(L)>0); % get rid of the zillion blank lines
ix=find(startsWith(L,'load')); % beginning of useful code by inspection
L([1:11 ix:end]) % display without extended header comment section
And, indeed -- this shows us a possible problem with T being a hardcoded number
T=114;
We've no way to know what this actually represents (although guessing from prior code it's an expected number of elements in a time series). This number may match up with the other data sets; is it the right number for this data set as well?
Again, have to read the code and make sure it has been modified everywhere it needs to be for a different data set than what was developed for -- hardcoding in magic numbers like this is extremely risky.
If you're lucky, it may be as simple as fixing that...
3 Comments
"T = 114 was the value also for Italy and it worked, although it was wrong!"
As we've discussed before, just because you can get code to run doesn't necessarily mean it is the correct answer. Again, one has to look back into how the code implements what it is supposed to be calculating. And, it's always possible the original authors made mistakes, too, given the way the code is structured and the mixing of data into code as is the case here.
You'll just have to dig deeper, then, into what is the root cause here; see above for what the algebra will require for this calculation to yield a positive result.
Although the following code seems quite unusual and there doesn't seem to be any real justification for it -- note that it creates the t2 variable out of simply a fixed percentage of this length number; whatever it is--
T = 114;
ix1=round(T*0.15)
ix2=round(T*0.85)
t2v=ix1:ix2;
numel(t2v)
Where did the 0.15 and 0.85 fractions come from and what is the justification for those particular fractions? Seems quite arbitrary...
Hence, later on, the L1, L2 values are simply pulled from whatever the result of the prior estimation of the MLE values were over this vector of t -- there seems no real reason to expect the specific conditions will be met to require the sign of the returned computed difference be negative in order for chow to be positive.
You'll simply have to look into what was computed for the est struct and understand how it relates to prior cases to figure out why this case doesn't appear to match other ones -- but it certainly would appear just from the code to be pretty much pure luck on which way the end result sign would turn out to be; don't see any reason from a theoretical standpoint it could be shown to have to turn out that way.
If I had to guess, I'd say it's likely the MLE values aren't the expected sign.
result=est(1:1:T);
L=result.likelihood;
param=result.param;
t2v=round(T*0.15):1:round(T*0.85);
Andrews=zeros(length(t2v),1);
AP=0;
Nyb=0;
AP0=0;
for t2=t2v; t2,
vecdates=1:1:t2; result=est(vecdates); L1=result.likelihood;
vecdates=t2+1:1:T; result=est(vecdates); L2=result.likelihood;
chow=-2* (L-(L1*(t2/T)+(1-t2/T)*L2));
Andrews(t2-round(T*0.15)+1,1)=chow;
end
Jinv=1/length(t2v);
Is est an array or a function? My initial supposition was it was an array being read from the input file, but that doesn't seem to make sense that est(1:T) would then result in what appears must be a single value....so, then you'll need to dig into what it does.
But, let's look at the rest of the code (again, paste the code and format it so can be read conveniently. "Help us help you...").
In the end it's looking for the maximum of any of these values it calculated and, as noted, for any value multiplied by a negative number to be positive, that value must also be negative. Looking at the max here means that unless every calculated value is negative, the resultant SuplR value is going to be greater than or equal to zero.
L=readlines('test.m');
L=L(strlength(L)>0);
ix=find(contains(L,'for t2'));
L(ix-5:ix+15)
run('test')
We can't run your code for lack of input data.in order to see, however, and, more than likely would have to have the missing m file(s) as well.
Note there's also another hardcoded number that looks suspiciously like T again...knowing the dates of interest, this should be an easily calculatable number instead to avoid these issues.
Categories
Find more on Variables in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
