Model solution
Contents
Model solution#
In this exercise, we will go through the model solution.
Warning
Remember that ams.func
is the sintax that calls the file with the functions see before.
Packages#
First, calls numpy and my functions file so we can use it in our model solution.
import numpy as np
import ams_functions as ams
Parameters#
To start our solution, we need to determine the parameters. The names here can be a little bit confusion, but the age_max
indicate how many adding yeas the child can study. This means that we are simulated children that are 6 years old and can study until they are 16.
age_start = 6
age_max = 10
Creating matrices#
Before the model solution loop, we need to create empty matrices to store our results. I commented the code so you can see what is the matrix about.
# Probabilities of working
probs_w = np.empty((age_max,age_max))
probs_w[:] = np.nan
# Probabilities of school
probs_s = np.empty((age_max,age_max))
probs_s[:] = np.nan
# Threshold matrix
eps_t = np.empty((age_max,age_max))
eps_t[:] = np.nan
# Expected values matrix for L = 1
eve_w = np.empty((age_max,age_max))
eve_w[:] = np.nan
# Expected values matrix for L = 0
eve_s = np.empty((age_max,age_max))
eve_s[:] = np.nan
# Value function matrix
vf = np.empty((age_max,age_max))
vf[:] = np.nan
Model solution loop#
The model solution is solved backward and follows these steps:
Solve for the threshold that makes the value function for going to school bigger.
Use the result from 1 to obtain the probability of going to school.
Calculate the truncated distribution value for going to school
Calculate the EV (expected values in the period ahead)
Repeat until the end of the loop
In this case, for the last period we will use the terminal value defined as:
This gives the value of having an specific level of education in the last period
for age in range(age_max,-1,-1):
for edu in range(0,age_max):
if edu >= age:
pass
else:
# For the threshold we will use the terminal value to capture the value of an specifiy education level in terms of utility
if age == age_max: #last period situation
tv_s = ams.prob_progress(edu)*(ams.terminal_v(edu = edu+1)) + (1 - ams.prob_progress(edu))*(ams.terminal_v(edu = edu))
tv_w = ams.terminal_v(edu = edu)
eps_t[age-1,edu] = ams.value_function(age,edu,1, EV = tv_s) - ams.value_function(age,edu,0, EV = tv_w)
else: #Other periods will take the difference between value functions
tv_s = ams.prob_progress(edu)*(vf[age, edu+1]) + (1 - ams.prob_progress(edu))*(vf[age,edu])
tv_w = vf[age, edu]
eps_t[age-1,edu] = ams.value_function(age,edu,1, EV = tv_s) - ams.value_function(age,edu,0, EV = tv_w)
# Now we will calcute the probabilities of schooling and working based on the thresholds
probs_s[age-1,edu] = ams.logistic(eps = eps_t[age-1,edu])
probs_w[age-1,edu] = 1 - probs_s[age-1,edu]
# Using the thresholds we will calculate the EVs
if age == age_max:
u_t = ams.trunc_change(eps = eps_t[age-1,edu])
eve_s[age-1,edu] = ams.value_function(age,edu,1, EV = tv_s) + ams.trunc_school(ϵ_threshold = eps_t[age-1,edu], u_threshold = u_t)
eve_w[age-1,edu] = ams.value_function(age,edu,0, EV = tv_w)
else:
u_t = ams.trunc_change(eps = eps_t[age-1,edu])
eve_s[age-1,edu] = ams.value_function(age,edu,1, EV = tv_s) + ams.trunc_school(ϵ_threshold = eps_t[age-1,edu], u_threshold = u_t)
eve_w[age-1,edu] = ams.value_function(age,edu,0, EV = tv_w)
vf[age-1,edu] = probs_s[age-1,edu]*eve_s[age-1,edu] + probs_w[age-1,edu]*eve_w[age-1,edu]
Let’s some results. The probability of going to school matrix looks like this:
probs_s
array([[0.99999993, nan, nan, nan, nan,
nan, nan, nan, nan, nan],
[0.99999884, 0.99999981, nan, nan, nan,
nan, nan, nan, nan, nan],
[0.99997971, 0.99999692, 0.99999955, nan, nan,
nan, nan, nan, nan, nan],
[0.99962732, 0.99994308, 0.99999232, 0.99999898, nan,
nan, nan, nan, nan, nan],
[0.99586383, 0.99880728, 0.99984954, 0.99998211, 0.99999649,
nan, nan, nan, nan, nan],
[0.91581254, 0.98644322, 0.9962496 , 0.99962788, 0.9999358 ,
0.99998729, nan, nan, nan, nan],
[0.65754572, 0.76395663, 0.96031838, 0.98826779, 0.99857002,
0.99975688, 0.99994754, nan, nan, nan],
[0.33117145, 0.39518156, 0.55636764, 0.91813627, 0.94320286,
0.99389726, 0.9989383 , 0.99971736, nan, nan],
[0.13441563, 0.15878852, 0.11884147, 0.7362356 , 0.67868231,
0.76872596, 0.97157708, 0.99372576, 0.99781425, nan],
[0.05324222, 0.05684499, 0.03267233, 0.43056996, 0.28367671,
0.27255018, 0.50621805, 0.88952229, 0.94513499, 0.97851892]])
The EV matrix is:
vf
array([[154.2996938 , nan, nan, nan,
nan, nan, nan, nan,
nan, nan],
[147.40631343, 165.27423062, nan, nan,
nan, nan, nan, nan,
nan, nan],
[142.41339016, 160.14706687, 179.02221917, nan,
nan, nan, nan, nan,
nan, nan],
[139.62337579, 157.12161177, 175.85698874, 195.86431707,
nan, nan, nan, nan,
nan, nan],
[139.31015202, 156.54370737, 175.01510424, 194.87866649,
211.0367933 , nan, nan, nan,
nan, nan],
[141.18076829, 158.80179566, 176.88834836, 196.46108271,
212.47273586, 228.6402713 , nan, nan,
nan, nan],
[146.3594149 , 163.55689397, 182.03683495, 201.05416816,
216.74777283, 232.76714143, 248.15710453, nan,
nan, nan],
[153.3799503 , 171.87520601, 189.96149747, 209.48071548,
224.37003875, 240.03721847, 255.27493485, 268.60861758,
nan, nan],
[161.33329164, 181.63603512, 201.23849672, 220.60613972,
236.52013511, 251.09686622, 265.88547942, 279.05383333,
289.06793576, nan],
[169.66802391, 192.06435084, 213.670688 , 233.8340728 ,
251.59177304, 267.19573071, 280.76515423, 293.43892839,
303.18467998, 311.39752683]])
Using our function file, you can use the following code to run the loop we just saw
probs_w, probs_s, eps_t, eve_w, eve_s, vf = ams.solution()