<jupyter_start><jupyter_text>Inaugural Project > **Note the following:** > 1. This is an example of how to structure your **inaugural 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 `inauguralproject.py` file includes a function which can be used multiple times in this notebook.  Getting an insuranceWe will here work with a benchmark model of insurance policies, Mossin's model. We are dealing with an agent wanting to get an insurance. The agent holds assets $y$ and may suffer a monetary loss $x$. For now we assume that the loss is incurred with probability $p$. Utility of assets is given by $u(\cdot)$.If the agent is not insured, expected value is$$V_{0}=p u(y-x)+(1-p) u(y)$$The agent may however obtain an insurance contract to mitigate a loss. The contract consists of a coverage amount $q$ and a premium $\pi$. The coverage cannot exceed the loss so $q \in[0, x]$ and $\pi$ must be paid irrespective of whether a loss was incurred. To not go broke, the insurance company must require higher premium for better coverage. It therefore implements a premium policy$$\pi(p, q)=p q,$$An insured agent will thus have expected utility$$V(q ; \pi)=p u(y-x+q-\pi(p, q))+(1-p) u(y-\pi(p, q))$$To solve the model numerically, we let $u$ be given by the function$$u(z)=\frac{z^{1+\vartheta}}{1+\vartheta}$$And use the following parameterization.$$y=1, p=0.2, \vartheta=-2$$ Imports and set magics:<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
plt.style.use('seaborn-whitegrid')
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm # for colormaps
from scipy import optimize
from scipy.stats import beta # beta distribution
# autoreload modules when code is run. Otherwise, python will not see recent changes. 
%load_ext autoreload
%autoreload 2
%reload_ext autoreload

# Import your own code
from inauguralproject import *<jupyter_output><empty_output><jupyter_text>Question 1 **Explain how you solve the model** *1. i) Construct a function that takes $(x, y, p)$ as arguments and returns the agents' optimal insurance coverage:$$q^{\star}=\operatorname{argmax}_{q \in[0, x]} V(q ; \pi)$$ ii) Make a grid of $x$ in the range $[0.01,0.9]$.<jupyter_code>N = 100 # I do 10 trials 
x_vec = np.linspace(0.01,0.9,N) # Here we have 10 values between 0.01 and 0.9
q_vec = [] # I make an empty list that<jupyter_output><empty_output><jupyter_text>iii) For each $x$, calculate $q^{\star}$.<jupyter_code>for i in x_vec: # Loop that goes through the 100 values in x_vec
    quess = 0
    obj = lambda q: - expected_utility_insured(i, q)
    res = optimize.minimize_scalar(obj, quess, method='bounded', bounds=(1e-8,i)) # Optimize with bounds from 0 to the element of x_vec. 
    q_op = res.x # The result is stored in the variable q_op
    q_vec.append(q_op)<jupyter_output><empty_output><jupyter_text>iv) Plot the $x$ s and $q^{\star} s$ and write what the graph suggests about the general solution to the agent's problem? (You can verify the finding analytically if you will).<jupyter_code>fig = plt.figure() 
ax = fig.add_subplot(1,1,1) # an empty canvas is created 

ax.plot(x_vec,q_vec,lw=2,color='black',label='$q$')

ax.set_xlabel('x - monetary loss')
ax.set_ylabel('q* - optimal coverage')
ax.legend(loc='lower center',frameon=True);<jupyter_output><empty_output><jupyter_text>It is clear that there is a positive linear relationship between $x$ and $q*$. This is in line with the analytical where to optimal coverage is $q*=x$  Question 2 We would like to find the set of acceptable contracts from the agent's point of view. This would be all the contracts $(q, \pi)$ which ensure an expected value at least as good as if not having an insurance, $V(q ; \pi)=V_{0}$. i) Let the loss be $x=0.6$ and construct a grid of $q$ over $[0.01,0.6]$. ii) Think now of $\pi$ as a variabel that has to be solved for (instead of it being a function of $p$ and $q$ as before). Loop over each element in the grid of $q$ s and find the corresponding $\tilde{\pi}$ such that $V(q ; \tilde{\pi})=V_{0}$.<jupyter_code>N = 100 # I do 10 trials 
q_grid = np.linspace(0.01,0.6,N) # Here we have 100 values between 0.01 and 0.9
pi_vec = np.empty(N) # I make an empty list that 
pi_quess = 0.0001

for i, q in enumerate(q_grid): # Loop that goes through the 100 values in 
    res = optimize.root(indiff, pi_quess, args=(q))
    pi_vec[i] = res.x[0]
print(pi_vec) # check if output looks correct<jupyter_output>[0.00605971 0.00963485 0.01318277 0.01670332 0.02019636 0.02366175
 0.02709936 0.03050907 0.03389076 0.03724432 0.04056965 0.04386667
 0.04713527 0.05037538 0.05358693 0.05676985 0.05992409 0.0630496
 0.06614633 0.06921426 0.07225336 0.07526361 0.07824501 0.08119755
 0.08412125 0.08701611 0.08988216 0.09271943 0.09552797 0.09830782
 0.10105904 0.10378168 0.10647583 0.10914156 0.11177895 0.11438811
 0.11696913 0.11952212 0.12204719 0.12454447 0.12701409 0.12945617
 0.13187087 0.13425833 0.13661871 0.13895216 0.14125884 0.14353894
 0.14579262 0.14802007 0.15022147 0.15239701 0.15454689 0.1566713
 0.15877044 0.16084454 0.16289378 0.16491839 0.16691859 0.16889459
 0.17084661 0.17277488 0.17467962 0.17656107 0.17841946 0.18025501
 0.18206796 0.18385855 0.18562702 0.18737359 0.18909852 0.19080203
 0.19248436 0.19414577 0.19578648 0.19740674 0.19900679 0.20058686
 0.2021472  0.20368804 0.20520962 0.20671218 0.20819596 0.20966119
 0.21110811 0.21253694 0.21394793 0.2153413  0.21671728 0.218076[...]<jupyter_text>iii) Make a plot of the acceptable premiums $\tilde{\pi}$ (the $q$ grid on the x-axis). Add to the diagram the function $\pi(p, q)=p q$ from before, which gives the premiums that the insurance company need at least in order to break even at each $q$. Taken together, these two graphs map out the set of feasible premiums for a given $q$.Try to make the graph pretty and informative (as if it was going into a text book or paper). That is, use shading, labels, titles etc.<jupyter_code>fig = plt.figure() 
ax = fig.add_subplot(1,1,1) # an empty canvas is created 

ax.plot(q_grid, pi_vec,lw=2,color='red',label='acceptable premuims for the agent')
ax.plot(q_grid, q_grid*0.2 ,lw=2,color='black',label='break even premiums for the company')
ax.fill_between(q_grid, q_grid*0.2, pi_vec  , color='grey', label = 'acceptable premiums for both')
ax.set_xlabel('q - grid')
ax.set_ylabel('$\pi$')
ax.legend(loc='lower center',frameon=True);<jupyter_output><empty_output><jupyter_text>Question 3 3. We will now consider a modification to the setup where the loss $x$ is drawn from a beta distribution and the coverage $q$ is a fraction of $x$. Specifically, we'll set$$\begin{gathered}q=\gamma x, \quad \gamma \in[0,1] \\x \sim \operatorname{Beta}(\alpha, \beta), \quad \alpha=2, \beta=7\end{gathered}$$which means that the agent's value is written as$$V(\gamma, \pi)=\int_{0}^{1} u(y-(1-\gamma) x-\pi) f(x) d x$$where $f(\cdot)$ is the density of the Beta distribution. Note that the Beta distribution has support in $[0,1]$ which explains the limits in the integral.<jupyter_code>np.random.seed(1) # We set a seed so the results can be replicated by others running the code. 
gammas = [0.9, 0.45] #We set the parameter values for gamma an pi. We then plug it into the function V
pies = [0.2, 0.1]

print("Consumer 1 =",  V(gammas[0], pies[0]))
print("Consumer 2 =",  V(gammas[1], pies[1]))<jupyter_output>Consumer 1 = -1.2858791338460491
Consumer 2 = -1.297714256088189<jupyter_text>Question 4 We now consider the situation from the insurance company’s viewpoint. They know theparameters of the loss distribution. Their policy is to let customers set a coverage ratioγ and then return an offer on the premium π to maximize profits. i) Given a customer wanting $\gamma$ = 0.95, what is the profit maximizing premium, $\pi^*$?<jupyter_code>pi_guess_1 = [0] 
obj = lambda pi:  V(pi[0],0.95) - V(0,0)
res_4 = optimize.root(obj, pi_guess_1, method = "broyden1") # We optimize and find the root

pi_res = res_4.x[0]

print("The profit maximizing premium is given by pi* =", pi_res)<jupyter_output>The profit maximizing premium is given by pi* = 0.5027297800503328
{'fraction_non_alphanumeric': 0.11369294605809128, 'fraction_numerical': 0.11867219917012448, 'mean_word_length': 4.654155495978553, 'pattern_counts': {'":': 0, '<': 30, '<?xml version=': 0, '>': 36, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '8', 'dir': 'jupyter-structured-clean-dedup', 'id': '66429', 'license': 'MIT', 'n_tokens_mistral': 3547, 'n_tokens_neox': 2857, 'n_words': 1076, 'path': 'inauguralproject/inauguralproject.ipynb', 'repo_name': 'NumEconCopenhagen/projects-2022-teamp-up', 'size': '8435'}