Skip to content

Step 1 of Jason's doc #64

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
evan-magnusson merged 1 commit into from
Jun 22, 2015
Merged

Step 1 of Jason's doc #64

Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
33 changes: 17 additions & 16 deletions Python/DEVELOPMENT/SS.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -138,10 +138,9 @@ def Euler_equation_solver(guesses, r, w, T_H, factor, j, params, chi_b, chi_n, t
error1[mask4] += 1e14
return list(error1.flatten()) + list(error2.flatten())

def wrguess(X, bssmat, nssmat, j, params, chi_b, chi_n, tau_bq, rho, lambdas, weights, e):
def wrguess(X, bssmat, nssmat, params, chi_b, chi_n, tau_bq, rho, lambdas, weights, e):

J, S, T, beta, sigma, alpha, Z, delta, ltilde, nu, g_y, tau_payroll, retire, mean_income_data, a_tax_income, b_tax_income, c_tax_income, d_tax_income, h_wealth, p_wealth, m_wealth, b_ellipse, upsilon = params

w, r, T_H, factor = X

for j in xrange(J):
Expand All @@ -153,11 +152,11 @@ def wrguess(X, bssmat, nssmat, j, params, chi_b, chi_n, tau_bq, rho, lambdas, we
# print np.array(Euler_equation_solver(np.append(bssmat[:, j], nssmat[:, j]), r, w, T_H, factor, j, params, chi_b, chi_n, theta, tau_bq, rho, lambdas, e)).max()

#Calculate demand for C, I, and Y
BQ = (weights * rho*bssmat).sum(0)
BQ = (weights * rho.reshape(S,1)*bssmat).sum(0).reshape(1,J)
theta = tax.replacement_rate_vals(nssmat, w, factor, e, J, weights)
net_tax = tax.total_taxes(r, bssmat, w, e, nssmat, BQ, lambdas, factor, T_H, 'SS', False, params, theta, tau_bq)
b_s = np.append(np.zeros(J).reshape(1,J), bssmat[:-1])
cons = house.get_cons(r, b_s, w, e, n, BQ, lambdas, bssmat, params, net_tax)
net_tax = tax.total_taxes(r, bssmat, w, e, nssmat, BQ, lambdas, factor, T_H, 0, 'SS', False, params, theta, tau_bq)
b_s = np.array(list(np.zeros(J).reshape(1,J)) + list(bssmat[:-1]))
cons = house.get_cons(r, b_s, w, e, nssmat, BQ, lambdas.reshape(1,J), bssmat, params, net_tax)
C = (weights*cons).sum()
B = (weights*bssmat).sum()
I = delta*B
Expand All @@ -172,12 +171,14 @@ def wrguess(X, bssmat, nssmat, j, params, chi_b, chi_n, tau_bq, rho, lambdas, we
error2 = L_demand - firm.get_L(e, nssmat, weights)

#Get other 2 errors from the factor equation and government constraint
revenue = tax.total_taxes(r, bssmat, w, e, nssmat, BQ, lambdas, factor, 0, 'SS', False, params, theta, tau_bq)
error3 = tax.total_taxes(r, bssmat, w, e, nssmat, BQ, lambdas, factor, 0, 0, 'SS', False, params, theta, tau_bq)
error3 = T_H - (weights*revenue).sum()
mean_income_model = ((r * b_s + w * e * nssmat) * weights).sum()
error4 = factor - mean_income_data / mean_income_model

return [error1, error2, error3, error4]
error = [error1, error2, error3, error4]
print error
return error


def SS_solver(b_guess_init, n_guess_init, wguess, rguess, T_Hguess, factorguess, chi_n, chi_b, params, iterative_params, tau_bq, rho, lambdas, weights, e):
Expand All @@ -189,7 +190,7 @@ def SS_solver(b_guess_init, n_guess_init, wguess, rguess, T_Hguess, factorguess,
bssmat = b_guess_init
nssmat = n_guess_init

w, r, T_H, factor = opt.fsolve(wrguess, [w,r,T_H,factor], args=(bssmat, nssmat, j, params, chi_b, chi_n, tau_bq, rho, lambdas, weights, e))
w, r, T_H, factor = opt.fsolve(wrguess, [w,r,T_H,factor], args=(bssmat, nssmat, params, chi_b, chi_n, tau_bq, rho, lambdas, weights, e))

eul_errors = np.ones(J)
b_mat = np.zeros((S, J))
Expand Down Expand Up @@ -288,15 +289,15 @@ def function_to_minimize(chi_params_scalars, chi_params_init, params, weights_SS
, 28.90029708 , 29.83586775 , 30.87563699 , 31.91207845 , 33.07449767
, 34.27919965 , 35.57195873 , 36.95045988 , 38.62308152])
chi_params[J:] = chi_n_guess
chi_params = list(chi_params)
# chi_params = list(chi_params)
# First run SS simulation with guesses at initial values for b, n, w, r, etc
b_guess = np.ones((S, J)).flatten() * .01
n_guess = np.ones((S, J)).flatten() * .5 * ltilde
wguess = 1.2
rguess = .06
T_Hguess = 0
factorguess = 100000
solutions = SS_solver(b_guess.reshape(S, J), n_guess.reshape(S, J), wguess, rguess, T_Hguess, factorguess, chi_params[J:], chi_params[:J], parameters, iterative_params, tau_bq, rho, lambdas, omega_SS, e)
solutions = SS_solver(b_guess.reshape(S, J), n_guess.reshape(S, J), wguess*.9, rguess*.9, T_Hguess*.9, factorguess*.9, chi_params[J:], chi_params[:J], parameters, iterative_params, tau_bq, rho, lambdas, omega_SS, e)
variables = ['solutions', 'chi_params']
dictionary = {}
for key in variables:
Expand All @@ -309,12 +310,12 @@ def function_to_minimize(chi_params_scalars, chi_params_init, params, weights_SS
# will be multiplied by the chi initial guesses inside the function. Otherwise, if chi^b_j=1e5 for some j, and the
# minimizer peturbs that value by 1e-8, the % difference will be extremely small, outside of the tolerance of the
# minimizer, and it will not change that parameter.
chi_params_scalars = np.ones(S+J)
chi_params_scalars = opt.minimize(function_to_minimize_X, chi_params_scalars, method='TNC', tol=1e-14, bounds=bnds, options={'maxiter': 1}).x
# chi_params_scalars = np.ones(S+J)
# chi_params_scalars = opt.minimize(function_to_minimize_X, chi_params_scalars, method='TNC', tol=1e-14, bounds=bnds, options={'maxiter': 1}).x
# chi_params_scalars = opt.minimize(function_to_minimize_X, chi_params_scalars, method='TNC', tol=1e-14, bounds=bnds).x
chi_params *= chi_params_scalars
print 'The final scaling params', chi_params_scalars
print 'The final bequest parameter values:', chi_params
# chi_params *= chi_params_scalars
# print 'The final scaling params', chi_params_scalars
# print 'The final bequest parameter values:', chi_params

solutions_dict = pickle.load(open("OUTPUT/Saved_moments/SS_init_solutions.pkl", "r"))
solutions = solutions_dict['solutions']
Expand Down