|
| 1 | +######## Solve |
| 2 | +import pandas |
| 3 | +import numpy |
| 4 | +from collections import OrderedDict |
| 5 | +# from solution import * |
| 6 | +from calibrations import * |
| 7 | +from time_consistent import solve as solve_time_consistent |
| 8 | +from time_consistent import simulate |
| 9 | +from bttt.trees import DeterministicTree, get_ts |
| 10 | +from matplotlib import pyplot as plt |
| 11 | +from bttt.trees import DeathTree |
| 12 | +from bttt.model import import_tree_model |
| 13 | +from matplotlib import pyplot as plt |
| 14 | +from numpy import * |
| 15 | + |
| 16 | + |
| 17 | + |
| 18 | +N = 25 |
| 19 | +T = 25 |
| 20 | + |
| 21 | + |
| 22 | +calib0 = calib.copy() |
| 23 | +beta = calib0['beta'] |
| 24 | +a = calib0['a'] |
| 25 | +p = calib0['p'] |
| 26 | +zbar = calib0['zbar'] |
| 27 | +max_R=5 |
| 28 | + |
| 29 | + |
| 30 | +tree = DeterministicTree(N) |
| 31 | +for s in tree.nodes: |
| 32 | + tree.values[s] = zbar |
| 33 | + model = import_tree_model('models.yaml', key='stochastic', tree=tree) |
| 34 | + |
| 35 | +def solve_it(**cc): |
| 36 | + tree = DeterministicTree(N) |
| 37 | + for s in tree.nodes: |
| 38 | + tree.values[s] = zbar |
| 39 | + model = import_tree_model('models.yaml', key='optimal', tree=tree) |
| 40 | + model.calibration.update(cc) |
| 41 | + sol = model.solve(verbose=True) |
| 42 | + df = numpy.concatenate( [get_ts(tree, sol, varname)[:,None] for varname in ['e','f', 'Gamma']], axis=1 ) |
| 43 | + df = pandas.DataFrame(df, columns=['e','f','Gamma']) |
| 44 | + return df |
| 45 | + |
| 46 | +from collections import OrderedDict |
| 47 | + |
| 48 | + |
| 49 | +Rvec0 = linspace(0.01, 3.0, 20) |
| 50 | + |
| 51 | + |
| 52 | + |
| 53 | +pvec = [0.8, 0.85, 0.9, 0.95, 1.0] |
| 54 | + |
| 55 | +# Unconstrained solutions |
| 56 | + |
| 57 | +unconstrained_sols = OrderedDict() |
| 58 | +for p in pvec: |
| 59 | + dfs = [solve_it(Rbar=i, min_f=0, p=p) for i in Rvec0] |
| 60 | + unconstrained_sols[p] = dfs |
| 61 | + |
| 62 | +all_gammas = OrderedDict() |
| 63 | +for p in pvec: |
| 64 | + dfs = unconstrained_sols[p] |
| 65 | + max_gammas = numpy.array([dfs[i]['Gamma'].max() for i,r in enumerate(Rvec0)]) |
| 66 | + all_gammas[p] = max_gammas |
| 67 | + |
| 68 | +for p in pvec: |
| 69 | + plt.plot(Rvec0, all_gammas[p]) |
| 70 | + |
| 71 | +alpha = 0.5 |
| 72 | + |
| 73 | +Rlimit = OrderedDict() |
| 74 | + |
| 75 | +for p in pvec: |
| 76 | + dfs = unconstrained_sols[p] |
| 77 | + max_gammas = numpy.array([dfs[i]['Gamma'].max() for i,r in enumerate(Rvec0)]) |
| 78 | + j = numpy.where(max_gammas<alpha)[0][0] |
| 79 | + a0 = max_gammas[j-1] |
| 80 | + a1 = max_gammas[j] |
| 81 | + r0 = Rvec0[j-1] |
| 82 | + r1 = Rvec0[j] |
| 83 | + rmax = r0 + (r1-r0)*(alpha-a0)/(a1-a0) |
| 84 | + Rlimit[p] = rmax |
| 85 | + |
| 86 | +# Constrained solutions |
| 87 | + |
| 88 | +R0 = 0.8 |
| 89 | + |
| 90 | +constrained_solutions = OrderedDict() |
| 91 | + |
| 92 | +for p in pvec: |
| 93 | + Rm = Rlimit[p] |
| 94 | + R1 = min([R0, Rm]) |
| 95 | + if Rm>0: |
| 96 | + df = solve_it(p=p, Rbar=R1) |
| 97 | + elif Rm<=0: |
| 98 | + df = solve_it(p=p, Rbar=0.001) |
| 99 | + df['R'] = R0 - df['f'].cumsum().shift() |
| 100 | + df['R'][0] = R0 |
| 101 | + constrained_solutions[p] = df |
| 102 | + |
| 103 | +d = {'dfs':constrained_solutions} |
| 104 | + |
| 105 | +import pickle |
| 106 | +with open("precomputed_alpha_p.pickle",'wb') as f: |
| 107 | + pickle.dump(d,f) |
| 108 | + |
| 109 | +### |
| 110 | +## |
| 111 | +### |
| 112 | + |
| 113 | + |
| 114 | +Rmax = 0.998 |
| 115 | +nRvec0 = numpy.array([0.01, 0.5, 0.9, Rmax, 2.0, 3.0]) |
| 116 | +ndfs = [solve_it(Rbar=i, min_f=0) for i in nRvec0] |
| 117 | +for i in [4,5]: |
| 118 | + df = ndfs[3].copy() |
| 119 | + ndfs[i] = df |
| 120 | +for i,df in enumerate(ndfs): |
| 121 | + df['R'] = nRvec0[i] - df['f'].cumsum() |
| 122 | +# |
| 123 | +# |
| 124 | +# # In[43]: |
| 125 | +# |
| 126 | +# |
| 127 | +# def plot_dfs(dfs, labels): |
| 128 | +# attributes = ['b', 'g', 'r'] |
| 129 | +# if not isinstance(dfs, list): |
| 130 | +# dfs = [dfs] |
| 131 | +# fig = plt.figure(figsize=(15,10)) |
| 132 | +# # plt.clear() |
| 133 | +# plt.subplot(131) |
| 134 | +# plt.plot(dfs[0].index,dfs[0].index*0+0.5,color='black', linestyle='--') |
| 135 | +# for i,df in enumerate(dfs): |
| 136 | +# plt.plot(df["Gamma"]) |
| 137 | +# |
| 138 | +# plt.grid() |
| 139 | +# plt.title("Marginal Value of Intervention") |
| 140 | +# # plt.figure() |
| 141 | +# plt.subplot(132) |
| 142 | +# for i,df in enumerate(dfs): |
| 143 | +# plt.plot(df["e"], label=labels[i]) |
| 144 | +# plt.legend(loc='lower right') |
| 145 | +# plt.grid() |
| 146 | +# plt.title("Exchange Rate") |
| 147 | +# plt.subplot(133) |
| 148 | +# for i,df in enumerate(dfs): |
| 149 | +# plt.plot(df["R"], label=labels[i]) |
| 150 | +# plt.grid() |
| 151 | +# plt.title("Reserves") |
| 152 | +# return fig |
| 153 | +# |
| 154 | +# |
| 155 | +# # In[44]: |
| 156 | + |
| 157 | + |
| 158 | +import pickle |
| 159 | +with open("precomputed_option_value.pickle","wb") as f: |
| 160 | + pickle.dump({"commitment":ndfs},f) |
| 161 | + |
| 162 | +# |
| 163 | +# # In[82]: |
| 164 | +# |
| 165 | +# |
| 166 | +# f = plot_dfs(ndfs, labels=['$R_0={}$'.format(i) for i in Rvec0]) |
| 167 | +# plt.savefig("optimal_choice_noconstraint.png", bbox_inches="tight") |
| 168 | +# plt.savefig("optimal_choice_noconstraint.pdf", bbox_inches="tight") |
| 169 | +# f |
| 170 | +# |
| 171 | +# |
| 172 | +# # In[11]: |
| 173 | +# |
| 174 | +# |
| 175 | +# Rvec0 = linspace(0.3, 3.0, 10) |
| 176 | +# dfs = [solve_it(Rbar=i) for i in Rvec0] |
| 177 | +# |
| 178 | +# |
| 179 | +# # In[12]: |
| 180 | +# |
| 181 | +# |
| 182 | +# f = plot_dfs(dfs, labels=['$R_0={}$'.format(i) for i in Rvec0]) |
| 183 | +# plt.savefig("optimal_choice_constraint.png", bbox_inches="tight") |
| 184 | +# plt.savefig("optimal_choice_constraint.pdf", bbox_inches="tight") |
| 185 | +# f |
0 commit comments