@@ -40,7 +40,7 @@ def get_L(n, p, method):
4040 L_presum = p .e [- 1 , :, :] * p .omega_SS * n
4141 L = L_presum .sum ()
4242 elif method == "TPI" :
43- L_presum = ( n * (p .e * p .omega [: p .T , :, :]) )
43+ L_presum = n * (p .e * p .omega [: p .T , :, :])
4444 L = L_presum .sum (1 ).sum (1 )
4545 return L
4646
@@ -67,28 +67,28 @@ def get_I(b_splus1, K_p1, K, p, method):
6767
6868 """
6969 if method == "SS" :
70- omega_extended = np .append (p .omega_SS [1 :, :], [0.0 ], axis = 0 )
71- imm_extended = np .append (p .imm_rates [- 1 , 1 :, :], [0.0 ], axis = 0 )
72- part2 = (
73- (
74- b_splus1 * omega_extended * imm_extended
75- ).sum ()
76- ) / (1 + p .g_n_ss )
70+ omega_extended = np .append (
71+ p .omega_SS [1 :, :], np .zeros ((1 , p .J )), axis = 0
72+ )
73+ imm_extended = np .append (
74+ p .imm_rates [- 1 , 1 :, :], np .zeros ((1 , p .J )), axis = 0
75+ )
76+ part2 = ((b_splus1 * omega_extended * imm_extended ).sum ()) / (
77+ 1 + p .g_n_ss
78+ )
7779 aggI = (1 + p .g_n_ss ) * np .exp (p .g_y ) * (K_p1 - part2 ) - (
7880 1.0 - p .delta
7981 ) * K
8082 elif method == "TPI" :
81- omega_shift = np .append (p .omega [: p .T , 1 :, :], np .zeros ((p .T , 1 , p .J )), axis = 1 )
83+ omega_shift = np .append (
84+ p .omega [: p .T , 1 :, :], np .zeros ((p .T , 1 , p .J )), axis = 1
85+ )
8286 imm_shift = np .append (
8387 p .imm_rates [: p .T , 1 :, :], np .zeros ((p .T , 1 , p .J )), axis = 1
8488 )
85- part2 = (
86- (
87- b_splus1 * imm_shift * omega_shift
88- )
89- .sum (1 )
90- .sum (1 )
91- ) / (1 + np .squeeze (np .hstack ((p .g_n [: p .T - 1 ], p .g_n_ss ))))
89+ part2 = ((b_splus1 * imm_shift * omega_shift ).sum (1 ).sum (1 )) / (
90+ 1 + np .squeeze (np .hstack ((p .g_n [: p .T - 1 ], p .g_n_ss )))
91+ )
9292 aggI = (
9393 1 + np .squeeze (np .hstack ((p .g_n [: p .T - 1 ], p .g_n_ss )))
9494 ) * np .exp (p .g_y ) * (K_p1 - part2 ) - (1.0 - p .delta ) * K
@@ -123,21 +123,31 @@ def get_B(b, p, method, preTP):
123123 if method == "SS" :
124124 if preTP :
125125 part1 = b * p .omega_S_preTP
126- omega_extended = np .append (p .omega_S_preTP [1 :, :], [0.0 ], axis = 0 )
127- imm_extended = np .append (p .imm_rates_preTP [1 :, :], [0.0 ], axis = 0 )
126+ omega_extended = np .append (
127+ p .omega_S_preTP [1 :, :], np .zeros ((1 , p .J )), axis = 0
128+ )
129+ imm_extended = np .append (
130+ p .imm_rates_preTP [1 :, :], np .zeros ((1 , p .J )), axis = 0
131+ )
128132 pop_growth_rate = p .g_n_preTP
129133 else :
130134 part1 = b * p .omega_SS
131- omega_extended = np .append (p .omega_SS [1 :, :], [0.0 ], axis = 0 )
132- imm_extended = np .append (p .imm_rates [- 1 , 1 :, :], [0.0 ], axis = 0 )
135+ omega_extended = np .append (
136+ p .omega_SS [1 :, :], np .zeros ((1 , p .J )), axis = 0
137+ )
138+ imm_extended = np .append (
139+ p .imm_rates [- 1 , 1 :, :], np .zeros ((1 , p .J )), axis = 0
140+ )
133141 pop_growth_rate = p .g_n_ss
134142 part2 = b * omega_extended * imm_extended
135143 B_presum = part1 + part2
136144 B = B_presum .sum ()
137145 B /= 1.0 + pop_growth_rate
138146 elif method == "TPI" :
139147 part1 = b * p .omega [: p .T , :, :]
140- omega_shift = np .append (p .omega [: p .T , 1 :, :], np .zeros ((p .T , 1 , p .J )), axis = 1 )
148+ omega_shift = np .append (
149+ p .omega [: p .T , 1 :, :], np .zeros ((p .T , 1 , p .J )), axis = 1
150+ )
141151 imm_shift = np .append (
142152 p .imm_rates [: p .T , 1 :, :], np .zeros ((p .T , 1 , p .J )), axis = 1
143153 )
@@ -183,21 +193,23 @@ def get_BQ(r, b_splus1, j, p, method, preTP):
183193 pop_growth_rate = p .g_n_ss
184194 rho = p .rho [- 1 , :]
185195 if j is not None :
186- BQ_presum = omega * rho * b_splus1
196+ BQ_presum = omega [:, j ] * rho [:, j ] * b_splus1
187197 else :
188198 BQ_presum = omega * rho * b_splus1
189199 BQ = BQ_presum .sum (0 )
190200 BQ *= (1.0 + r ) / (1.0 + pop_growth_rate )
191201 elif method == "TPI" :
192202 pop = np .append (
193- p .omega_S_preTP .reshape (1 , p .S ), p .omega [: p .T - 1 , :], axis = 0
203+ p .omega_S_preTP .reshape (1 , p .S , p .J ),
204+ p .omega [: p .T - 1 , :, :],
205+ axis = 0 ,
194206 )
195207 rho = np .append (
196- p .rho_preTP .reshape (1 , p .S ), p .rho [: p .T - 1 , :], axis = 0
208+ p .rho_preTP .reshape (1 , p .S , p . J ), p .rho [: p .T - 1 , : , :], axis = 0
197209 )
198210
199211 if j is not None :
200- BQ_presum = b_splus1 * pop [:, j ] * rho [:, j ]
212+ BQ_presum = b_splus1 * pop [:, :, j ] * rho [:, :, j ]
201213 BQ = BQ_presum .sum (1 )
202214 BQ *= (1.0 + r ) / (1.0 + np .append (p .g_n_preTP , p .g_n [: p .T - 1 ]))
203215 else :
@@ -289,13 +301,9 @@ def get_C(c, p, method):
289301 """
290302
291303 if method == "SS" :
292- aggC = (
293- (c * p .omega_SS ).sum (- 1 ).sum (- 1 )
294- )
304+ aggC = (c * p .omega_SS ).sum (- 1 ).sum (- 1 )
295305 elif method == "TPI" :
296- aggC = (
297- (c * p .omega [: p .T , :, :]).sum (- 1 ).sum (- 1 )
298- )
306+ aggC = (c * p .omega [: p .T , :, :]).sum (- 1 ).sum (- 1 )
299307 return aggC
300308
301309
0 commit comments