<jupyter_start><jupyter_text>YOUR PROJECT TITLE > **Note the following:** > 1. This is *not* meant to be an example of an actual **model analysis project**, just an example of how to structure such a project.> 1. Remember the general advice on structuring and commenting your code from [lecture 5](https://numeconcopenhagen.netlify.com/lectures/Workflow_and_debugging).> 1. Remember this [guide](https://www.markdownguide.org/basic-syntax/) on markdown and (a bit of) latex.> 1. Turn on automatic numbering by clicking on the small icon on top of the table of contents in the left sidebar.> 1. The `modelproject.py` file includes a function which could be used multiple times in this notebook. Imports and set magics:<jupyter_code>#import numpy as np
#from scipy import optimize
#import sympy as sm

# autoreload modules when code is run
#%load_ext autoreload
#%autoreload 2

# local modules
#import modelproject<jupyter_output><empty_output><jupyter_text>Model description **Write out the model in equations here.** Make sure you explain well the purpose of the model and comment so that other students who may not have seen it before can follow.    Analytical solution If your model allows for an analytical solution, you should provide here.You may use Sympy for this. Then you can characterize the solution as a function of a parameter of the model.To characterize the solution, first derive a steady state equation as a function of a parameter using Sympy.solve and then turn it into a python function by Sympy.lambdify. See the lecture notes for details.   Numerical solution You can always solve a model numerically. Define first the set of parameters you need. Then choose one of the optimization algorithms that we have gone through in the lectures based on what you think is most fitting for your model.Are there any problems with convergence? Does the model converge for all starting values? Make a lot of testing to figure these things out.   Further analysis Make detailed vizualizations of how your model changes with parameter values. Try to make an extension of the model.   Conclusion Add concise conclusion.  Writing latex using:https://personal.math.ubc.ca/~pwalls/math-python/jupyter/latex/ Brug scipy linalg til at løse ax=b for x.sympy kan bruges med symboler og beregne som man normalt gør. Dette skriver resultaterne som funktioner med tegn.   Modelproject  Importing<jupyter_code>import modelproject
import sympy as sm
import numpy as np
import ipywidgets as widgets
import matplotlib.pyplot as plt

%load_ext autoreload
%autoreload 2<jupyter_output><empty_output><jupyter_text>Describe the model and the intention **The model**In this project, we examine the OLG-model. It is based on an idea, that agents live for two periods, and maximizes their utility in these.Utility is given by: $$U_t=ln(c_{1t})+\beta ln(c_{2t+1})$$ where: $$\beta>0$$Consumption is given by:$$c_{1t}=w_t-s_t$$$$c_{2t+1}=R_{t+1} s_t$$Assuming a cobb douglas production function:$$Y_t=K_t^\alpha L_t^{1-\alpha}$$where: $$0<\alpha<1$$**Variables and parameters**In the model, small letters are pr. worker terms, and capital letters are agregate.U is utilityC is consumption. 1t is for young individuals in period t, while 2t+1 is for the old generation in periode t+1K is capitalS is savingW is wageR is interest ratebeta is the time preferencesalpha is the share of capital in production **Extension for model from PS5 OPG2:**$$Y_t^i=A(K_t^i)^\alpha (K_t N_t^i)^{1-\alpha}$$where: $$A>0$$ **Extension for model with PAYG og FF:****Pay as you go:**$$c_{1t}=w_t(1-\tau)-s_t$$$$c_{2t+1}=s_tR_{t+1}+{\tau}w_{t+1}$$**Fully funded:**$$c_{1t}=w_t(1-\tau)-s_t$$$$c_{2t+1}=(s_t+{\tau}w_t)R_{t+1}$$**Fully funded is always prefered, since it does not affect capital accumulation**  Analytical solution We start by defining the variables with sympy, in order to solve the model analytically:<jupyter_code># Aggregate variables
U_t = sm.symbols('U_t')
R_t1 = sm.symbols('R_t1')
Y_t = sm.symbols('Y_t')
K_t = sm.symbols('K_t')
L_t = sm.symbols('L_t')
N_t = sm.symbols('N_t')
A = sm.symbols('A')
pi_t = sm.symbols('Pi_t')

# pr. worker variables
r_t = sm.symbols('r_t')
r_t_star = sm.symbols('r_t^*')
r_star = sm.symbols('r^*')
s_t = sm.symbols('s_t')
s_t_star = sm.symbols('s_t^*')
s_star = sm.symbols('s^*')
w_t = sm.symbols('w_t')
w_star = sm.symbols('w^*')
k_t = sm.symbols('k_t')
k_ss = sm.symbols('k_ss')
k_star = sm.symbols('k^*')
c_1t = sm.symbols('c_1t')
c_2t1 = sm.symbols('c_2t+1')
i = sm.symbols('i')
n = sm.symbols('n')

# Parameters
alpha = sm.symbols('alpha')
beta = sm.symbols('beta')<jupyter_output><empty_output><jupyter_text>The firms profit is given by their production minus the cost of inputs:<jupyter_code># The firm maximizes profit
profit = K_t**alpha * L_t**(1-alpha) - r_t * K_t - w_t*L_t
sm.Eq(pi_t,profit)<jupyter_output><empty_output><jupyter_text>Now the profitmaximizing rent is found by finding the first order condition and then solve:<jupyter_code># The profitmaximizing rent is found
# FOC wrt. capital
profit_diff_K = sm.diff(profit,K_t)

# Solving for rental rent
r_solution = sm.solve(profit_diff_K,r_t)[0]
#r_solution.subs(K_t**(alpha-1)*L_t**(1-alpha),k_t**(alpha-1))
sm.Eq(r_t,r_solution)
# The rent rewritten
r_rewritten = alpha*k_t**(alpha-1)

# Saving the rent as a python function:
r_ss = sm.lambdify((alpha,k_t),r_rewritten)<jupyter_output><empty_output><jupyter_text>Now the profitmaximizing wage is found in the same way as the rent:<jupyter_code># The profitmaximizing wage is found
# FOC wrt. labor
profit_diff_L = sm.diff(profit,L_t)

# Solving for wage
w_solution = sm.solve(profit_diff_L,w_t)[0]
sm.Eq(w_t,w_solution)
# The wage rewritten
w_rewritten = (1-alpha)*k_t**alpha

# Saving the wage as a python function:
w_ss = sm.lambdify((alpha,k_t),w_rewritten)<jupyter_output><empty_output><jupyter_text>Given the wage and rent, the agents maximizes utility in the two periods:<jupyter_code># The utilitymaximizing saving is found
# The budgetconstraints is given by
budgetcons_1 = w_t - s_t
budgetcons_2 = s_t*R_t1

# The budget constraints are inserted into the utility
utility = sm.ln(budgetcons_1) + beta*sm.ln(budgetcons_2)
utility

# FOC wrt. s
utility_diff_s = sm.diff(utility,s_t)

# Solving for s
s_solution = sm.solve(utility_diff_s,s_t)[0]
sm.Eq(s_t,s_solution)
# The saving rewritten
s_rewritten = (beta/(beta+1))*(1-alpha)*k_t**alpha

# Saving the savings as a python function:
s_ss = sm.lambdify((beta,alpha,k_t),s_rewritten)<jupyter_output><empty_output><jupyter_text>We use the savings to find the staedy state capital. We have that:$$K_{t+1}=S_t$$and$$L_{t+1}={L_t}(1+n)$$<jupyter_code># The steady state capital is found
# In steady state we have k_t+1=s_t:
savings = ((beta)/((1+beta)*(1+n)))*(k_ss**alpha)*(1-alpha)-k_ss

k_solution = sm.solve(savings,k_ss)[0]
k_solution
sm.Eq(k_ss,k_solution)
# The capital rewritten:
k_rewritten = ((beta*(1-alpha))/((1+beta)*(1+n)))**(1/(1-alpha))

# Saving the capital as a python function:
k_ss = sm.lambdify((beta, alpha, n),k_rewritten)
# We have the following results:
display(sm.Eq(r_star,r_rewritten),
        sm.Eq(w_star,w_rewritten),
        sm.Eq(s_star,s_rewritten),
        sm.Eq(k_star,k_rewritten))<jupyter_output><empty_output><jupyter_text>Numerical solution **Solving nummerically**Now we solve the model nummerically. We assume that the agent place higher utility on consumption in the first period than on consumption in the second period. Therefore beta is below 1.<jupyter_code># The capital 
beta_estimat = 2/3
alpha_estimat = 1/3
n_estimat = 0.1
ss_capital = k_ss(beta_estimat,alpha_estimat,n_estimat)
ss_saving = s_ss(beta_estimat,alpha_estimat,ss_capital)
ss_rent = r_ss(alpha_estimat, ss_capital)
ss_wage = w_ss(alpha_estimat, ss_capital)

print("Rent in steady state is: ", ss_rent)
print("Wage in steady state is: ", ss_wage)
print("Saving in steady state is: ", ss_saving)
print("Capital in steady state is: ", ss_capital)<jupyter_output>Rent in steady state is:  1.3749999999999996
Wage in steady state is:  0.3282439759448874
Saving in steady state is:  0.13129759037795496
Capital in steady state is:  0.1193614457981409<jupyter_text>Monte carlo simulation:**Simulating utility**We start by simulating the utility as a function of beta, which we assume to be uniformly destributed.<jupyter_code>print(modelproject.MonteCarlo(alpha_estimat, n_estimat))<jupyter_output>0.1200501800356351<jupyter_text>Interactive graph<jupyter_code>widgets.interact(modelproject.fig, 
    beta = widgets.FloatSlider(description='beta', min=0, max=1, step=0.01, value=0.5),
    alpha = widgets.fixed(1/3)
);<jupyter_output><empty_output><jupyter_text>Extension**Pay as you go:**$$c_{1t}=w_t(1-\tau)-s_t$$$$c_{2t+1}=s_tR_{t+1}+{\tau}w_{t+1}$$  Analytical<jupyter_code>w_t1 = sm.symbols('w_t1')
k_ss_PAYG = sm.symbols('k_ss(PAYG)')
tau = sm.symbols('tau')
# The utilitymaximizing saving is found
# The budgetconstraints is given by
budgetcons_1_PAYG = w_t*(1-tau) - s_t
budgetcons_2_PAYG = s_t*R_t1+tau*w_t1

# The budget constraints are inserted into the utility
utility_PAYG = sm.ln(budgetcons_1_PAYG) + beta*sm.ln(budgetcons_2_PAYG)
utility_PAYG

# FOC wrt. s
utility_diff_s_PAYG = sm.diff(utility_PAYG,s_t)

# Solving for s
s_solution_PAYG = sm.solve(utility_diff_s_PAYG,s_t)[0]
sm.Eq(s_t,s_solution_PAYG)
# The saving rewritten
s_rewritten_PAYG = (((R_t1*beta*w_t)*(1-tau))-tau*w_t1)/(R_t1*(beta+1))

# Saving the savings as a python function:
s_ss_PAYG = sm.lambdify((beta,alpha,tau,R_t1,w_t,w_t1),s_rewritten_PAYG)
s_rewritten_PAYG
# The steady state capital is found
# In steady state we have k_t+1=s_t:
savings_PAYG = (((R_t1*beta*(1-alpha)*(k_ss_PAYG**alpha))*(1-tau))-tau*((1-alpha)*(k_ss_PAYG**alpha)))/(R_t1*(beta+1))-k_ss_PAYG

k_solution_PAYG = sm.solve(savings_PAYG,k_ss_PAYG)[0]
k_solution_PAYG
sm.Eq(k_ss_PAYG,k_solution_PAYG)
# The capital rewritten:
k_rewritten_PAYG = ((R_t1*beta*(alpha-1)*(tau-1)+alpha*tau-tau)/(R_t1*(beta+1)))**(1/(1-alpha))

# Saving the capital as a python function:
k_ss_PAYG = sm.lambdify((beta, alpha, n, tau, R_t1),k_rewritten_PAYG)
# We have the following results:
display(sm.Eq(r_star,r_rewritten),
        sm.Eq(w_star,w_rewritten),
        sm.Eq(s_star,s_rewritten_PAYG),
        sm.Eq(k_star,k_rewritten_PAYG))<jupyter_output><empty_output><jupyter_text>Numerical<jupyter_code># The capital 
beta_estimat = 2/3
alpha_estimat = 1/3
n_estimat = 0.1
tau_estimat = 0.2
R_t1_estimat = 1.2

ss_capital_PAYG = k_ss_PAYG(beta_estimat,alpha_estimat,n_estimat, tau_estimat, R_t1_estimat)
ss_rent_PAYG = r_ss(alpha_estimat, ss_capital_PAYG)
ss_wage_PAYG = w_ss(alpha_estimat, ss_capital_PAYG)
ss_saving_PAYG = s_ss_PAYG(beta_estimat,alpha_estimat,tau_estimat,R_t1_estimat,ss_wage_PAYG,ss_wage_PAYG)

print("Rent in steady state is with PAYG: ", ss_rent_PAYG)
print("Wage in steady state is with PAYG: ", ss_wage_PAYG)
print("Saving in steady state is with PAYG: ", ss_saving_PAYG)
print("Capital in steady state is with PAYG: ", ss_capital_PAYG)
s_rewritten_PAYG = ((((1.2)*(2/3)*w_t)*(1-(0.2)))-(0.2)*w_t1)/((1.2)*((2/3)+1))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12179082965498604, 'fraction_numerical': 0.021980001801639493, 'mean_word_length': 6.1395498392282954, 'pattern_counts': {'":': 0, '<': 53, '<?xml version=': 0, '>': 59, 'https://': 3, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 1, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '56012', 'license': 'MIT', 'n_tokens_mistral': 4404, 'n_tokens_neox': 4213, 'n_words': 1101, 'path': 'modelproject/modelproject.ipynb', 'repo_name': 'NumEconCopenhagen/projects-2022-jonathan-og-anna', 'size': '11101'}