求MS-MIDAS模型GAUSS代码

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

有人有MS-MIDAS模型的GAUSS代码吗？我这里有一部分，但我总感觉不全。

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：MIDAS GAUSS Aus Das MID

1. new; cls;
   
2. library maxlikmt;
   
3. 
   
4. /*
   
5.     MS-MIDAS-AR (2-state) template in GAUSS
   
6.     ---------------------------------------
   
7.     Model:
   
8.       y_t = c_{S_t} + phi * y_{t-1} + b_{S_t} * (X_t * w(theta)) + e_t
   
9.       e_t | S_t=s ~ N(0, sigma_s^2)
   
10. 
    
11.     Regime switching:
    
12.       P = [ p11      1-p11
    
13.             1-p22    p22  ]
    
14. 
    
15.     Assumptions:
    
16.       1) y      : T x 1 low-frequency dependent variable
    
17.       2) ylag   : T x 1 lagged y
    
18.       3) Xlags  : T x K matrix of aligned high-frequency lags
    
19.                   Each row already contains the K high-frequency lags
    
20.                   relevant for the low-frequency observation.
    
21.       4) Common MIDAS weights across regimes (easy to extend later)
    
22. */
    
23. 
    
24. 
    
25. /* ---------- helpers ---------- */
    
26. 
    
27. proc (1) = logistic(x);
    
28.     retp(1 ./ (1 + exp(-x)));
    
29. endp;
    
30. 
    
31. proc (1) = logdnorm(y, mu, sig);
    
32.     local z;
    
33.     z = (y - mu) ./ sig;
    
34.     retp(-0.5 * ln(2*pi) - ln(sig) - 0.5 * (z .* z));
    
35. endp;
    
36. 
    
37. proc (1) = expalmon(theta, K);
    
38.     local j, z, w;
    
39.     j = seqa(0, 1, K);              // 0,1,2,...,K-1
    
40.     z = theta[1] .* j + theta[2] .* (j .* j);
    
41.     z = z - maxc(z);                // stabilize exponentials
    
42.     w = exp(z);
    
43.     retp(w ./ sumc(w));
    
44. endp;
    
45. 
    
46. proc (1) = stationary_probs(p11, p22);
    
47.     local den, pi1;
    
48.     den = 2 - p11 - p22;
    
49.     pi1 = (1 - p22) / den;
    
50.     retp(pi1 | (1 - pi1));
    
51. endp;
    
52. 
    
53. 
    
54. /* ---------- log-likelihood ---------- */
    
55. 
    
56. proc lnlk_ms_midas(struct PV p, struct DS d, ind);
    
57.     local y, ylag, Xlags, T, K;
    
58.     local c1, c2, b1, b2, phi, ls1, ls2, th1, th2, lp11, lp22;
    
59.     local s1, s2, p11, p22, P, w, xw;
    
60.     local alpha, pred, ll, t, lpdf1, lpdf2, m, dens, ft, mu1, mu2;
    
61. 
    
62.     struct modelResults mm;
    
63. 
    
64.     y     = d[1].dataMatrix;
    
65.     ylag  = d[2].dataMatrix;
    
66.     Xlags = d[3].dataMatrix;
    
67. 
    
68.     T = rows(y);
    
69.     K = cols(Xlags);
    
70. 
    
71.     /* unpack parameters */
    
72.     c1   = pvUnpack(p, 1);
    
73.     c2   = pvUnpack(p, 2);
    
74.     b1   = pvUnpack(p, 3);
    
75.     b2   = pvUnpack(p, 4);
    
76.     phi  = pvUnpack(p, 5);
    
77.     ls1  = pvUnpack(p, 6);   // log sigma_1
    
78.     ls2  = pvUnpack(p, 7);   // log sigma_2
    
79.     th1  = pvUnpack(p, 8);   // exp-almon theta_1
    
80.     th2  = pvUnpack(p, 9);   // exp-almon theta_2
    
81.     lp11 = pvUnpack(p,10);   // logit(p11)
    
82.     lp22 = pvUnpack(p,11);   // logit(p22)
    
83. 
    
84.     s1  = exp(ls1);
    
85.     s2  = exp(ls2);
    
86.     p11 = logistic(lp11);
    
87.     p22 = logistic(lp22);
    
88. 
    
89.     w  = expalmon(th1 | th2, K);
    
90.     xw = Xlags * w;
    
91. 
    
92.     P = p11~(1-p11) |
    
93.         (1-p22)~p22;
    
94. 
    
95.     alpha = stationary_probs(p11, p22);
    
96.     ll = zeros(T, 1);
    
97. 
    
98.     for t (1, T, 1);
    
99. 
    
100.         if t > 1;
     
101.             pred = P' * alpha;
     
102.         else;
     
103.             pred = alpha;
     
104.         endif;
     
105. 
     
106.         mu1 = c1 + phi * ylag[t] + b1 * xw[t];
     
107.         mu2 = c2 + phi * ylag[t] + b2 * xw[t];
     
108. 
     
109.         lpdf1 = logdnorm(y[t], mu1, s1);
     
110.         lpdf2 = logdnorm(y[t], mu2, s2);
     
111. 
     
112.         /* log-sum-exp stabilization */
     
113.         m    = maxc(lpdf1 | lpdf2);
     
114.         dens = exp((lpdf1 | lpdf2) - m);
     
115.         ft   = sumc(pred .* dens);
     
116. 
     
117.         if ft <= 1e-300;
     
118.             if ind[1];
     
119.                 mm.function = -1e256 * ones(T,1);
     
120.             endif;
     
121.             retp(mm);
     
122.         endif;
     
123. 
     
124.         ll[t] = ln(ft) + m;
     
125.         alpha = (pred .* dens) / ft;
     
126.     endfor;
     
127. 
     
128.     if ind[1];
     
129.         mm.function = ll;   // observation-wise log-likelihood
     
130.     endif;
     
131. 
     
132.     retp(mm);
     
133. endp;
     
134. 
     
135. 
     
136. /* ---------- example estimation block ---------- */
     
137. 
     
138. /*
     
139.     Replace these with your actual aligned data:
     
140. 
     
141.     y      : T x 1
     
142.     ylag   : T x 1
     
143.     Xlags  : T x K
     
144. */
     
145. 
     
146. /* Example placeholders:
     
147. y     = loadd("y.dat");
     
148. ylag  = loadd("ylag.dat");
     
149. Xlags = loadd("Xlags.dat");
     
150. */
     
151. 
     
152. struct DS d0;
     
153. d0 = reshape(dsCreate, 3, 1);
     
154. d0[1].dataMatrix = y;
     
155. d0[2].dataMatrix = ylag;
     
156. d0[3].dataMatrix = Xlags;
     
157. 
     
158. struct PV p0;
     
159. p0 = pvPacki(pvCreate,  0.00, "c1",    1);
     
160. p0 = pvPacki(p0,        0.50, "c2",    2);
     
161. p0 = pvPacki(p0,        0.20, "b1",    3);
     
162. p0 = pvPacki(p0,        0.05, "b2",    4);
     
163. p0 = pvPacki(p0,        0.30, "phi",   5);
     
164. p0 = pvPacki(p0,       -1.00, "ls1",   6);
     
165. p0 = pvPacki(p0,       -0.50, "ls2",   7);
     
166. p0 = pvPacki(p0,        0.05, "th1",   8);
     
167. p0 = pvPacki(p0,       -0.01, "th2",   9);
     
168. p0 = pvPacki(p0,        2.00, "lp11", 10);
     
169. p0 = pvPacki(p0,        2.00, "lp22", 11);
     
170. 
     
171. struct maxlikmtControl c0;
     
172. c0 = maxlikmtControlCreate;
     
173. c0.printIters = 1;
     
174. c0.maxIters   = 5000;
     
175. c0.tol        = 1e-6;
     
176. /* optional: c0.algorithm = "BFGS"; */
     
177. 
     
178. struct maxlikmtResults out;
     
179. out = maxlikmt(&lnlk_ms_midas, p0, d0, c0);
     
180. 
     
181. print "retcode = " out.retcode;
     
182. print "logLik  = " out.fct;
     
183. print "estimated parameters:";
     
184. print out.par;

复制代码