Skip to content
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

HW1 submission #3

Open
wants to merge 2 commits into
base: master
Choose a base branch
from
Open

HW1 submission #3

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
ac567b0
partly completed HW
Mar 7, 2018
9e1f959
HW1 submission
Mar 7, 2018
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
270 changes: 270 additions & 0 deletions questions_AD.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,270 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Integration exercises\n",
"\n",
"1. Fackler and Miranda exercise 5.1: Demand for a commodity is given by $q(p)=2p^{-0.5}$. The price of the good falls from 4 to 1. Compute the change in consumer surplus **analytically** using calculus and **numerically** using three different methods: a Gauss-Legendre rule, Monte-Carlo, and Quasi Monte Carlo.\n",
" 1. Start out by plotting the demand function. In that plot, label the axis properly ($p$ on the $x$ axis), and add 2 horizontal lines for the equilibrium quantities at both prices $p=1,p=4$. This should help you visualize the consumer surplus, and it should guide you in the analytical solution. The simplest way to solve this first part would be to edit this `IJulia` notebook, just adding a code cell below (for the plot), and then one for the analytic solution. You can get and see how to use `IJulia` [by clicking on this link](https://github.com/JuliaLang/IJulia.jl).\n",
"\n",
" To implement the integration rules, write 3 functions, one for each sub question, and detailed below. Let each function take as argument the number of integration points. Use $n=10,15,20$ points. Each of those functions should produce the result (of course), as well as a plot where we can see integration nodes vs function value. We want to get an idea where each method places the points, and how this might influence the different results. Ideally, the result would tell us how far the corresponding method is away from your analytic solution. *Ideally*, your plot would show us all the results at once, so we can easily compare across methods. Keyword *subplots*.\n",
" 1. Gauss-legendre rule. Note that you have to change the function domain to $[-1,1]$ first. This is achieved with the following transformation\n",
"$$ \\int_a^b f(x)\\,dx = \\frac{b-a}{2} \\int_{-1}^1 f\\left(\\frac{b-a}{2}x + \\frac{a+b}{2}\\right)dx $$\n",
" 1. monte-carlo, again taking $n$ as an argument.\n",
" 1. pseudo monte-carlo. Use a Sobol Sequence.\n",
"2. Fackler and Miranda exercise 5.5: A government stabilizes the supply of a commodity at $S=2$, but allows the price to be determined on the market. Domestic and export demand for the commodity are given by $D = \\tilde{\\theta_1} p^{-1}, X = \\tilde{\\theta_2} p^{-0.5}$, and where $\\log \\tilde{\\theta_1}$ and $\\log \\tilde{\\theta_1}$ are both normally distributed with means 0, variances 0.02 and 0.01 respectively, and covariance 0.01.\n",
" 1. compute the expected price $E[p]$ and expected variance $V[p]$ using a kronecker product rule with gauss-hermite grids (10) in each dimension of the shock.\n",
" 1. perform the same computation with a monte carlo integration scheme.\n",
" 1. Bonus question: same thing with quasi-monte carlo."
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"#### Question 1. Part A."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<?xml version=\"1.0\" encoding=\"utf-8\"?>\n",
"<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"600\" height=\"400\" viewBox=\"0 0 600 400\">\n",
"<defs>\n",
" <clipPath id=\"clip1000\">\n",
" <rect x=\"0\" y=\"0\" width=\"600\" height=\"400\"/>\n",
" </clipPath>\n",
"</defs>\n",
"<polygon clip-path=\"url(#clip1000)\" points=\"\n",
"0,400 600,400 600,0 0,0 \n",
" \" fill=\"#ffffff\" fill-opacity=\"1\"/>\n",
"<defs>\n",
" <clipPath id=\"clip1001\">\n",
" <rect x=\"120\" y=\"0\" width=\"421\" height=\"400\"/>\n",
" </clipPath>\n",
"</defs>\n",
"<polygon clip-path=\"url(#clip1000)\" points=\"\n",
"48.2225,360.121 580.315,360.121 580.315,31.4961 48.2225,31.4961 \n",
" \" fill=\"#ffffff\" fill-opacity=\"1\"/>\n",
"<defs>\n",
" <clipPath id=\"clip1002\">\n",
" <rect x=\"48\" y=\"31\" width=\"533\" height=\"330\"/>\n",
" </clipPath>\n",
"</defs>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 106.934,360.121 106.934,31.4961 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 226.996,360.121 226.996,31.4961 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 347.058,360.121 347.058,31.4961 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 467.12,360.121 467.12,31.4961 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 48.2225,342.737 580.315,342.737 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 48.2225,256.178 580.315,256.178 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 48.2225,169.619 580.315,169.619 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#000000; stroke-width:0.5; stroke-opacity:0.1; fill:none\" points=\"\n",
" 48.2225,83.0602 580.315,83.0602 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,360.121 580.315,360.121 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,360.121 48.2225,31.4961 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 106.934,360.121 106.934,355.191 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 226.996,360.121 226.996,355.191 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 347.058,360.121 347.058,355.191 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 467.12,360.121 467.12,355.191 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,342.737 56.2039,342.737 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,256.178 56.2039,256.178 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,169.619 56.2039,169.619 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,83.0602 56.2039,83.0602 \n",
" \"/>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:middle;\" transform=\"rotate(0, 106.934, 373.921)\" x=\"106.934\" y=\"373.921\">1</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:middle;\" transform=\"rotate(0, 226.996, 373.921)\" x=\"226.996\" y=\"373.921\">2</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:middle;\" transform=\"rotate(0, 347.058, 373.921)\" x=\"347.058\" y=\"373.921\">3</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:middle;\" transform=\"rotate(0, 467.12, 373.921)\" x=\"467.12\" y=\"373.921\">4</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:end;\" transform=\"rotate(0, 42.2225, 347.237)\" x=\"42.2225\" y=\"347.237\">1.0</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:end;\" transform=\"rotate(0, 42.2225, 260.678)\" x=\"42.2225\" y=\"260.678\">1.5</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:end;\" transform=\"rotate(0, 42.2225, 174.119)\" x=\"42.2225\" y=\"174.119\">2.0</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:end;\" transform=\"rotate(0, 42.2225, 87.5602)\" x=\"42.2225\" y=\"87.5602\">2.5</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:21; text-anchor:middle;\" transform=\"rotate(0, 314.269, 18)\" x=\"314.269\" y=\"18\">Demand curve</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:16; text-anchor:middle;\" transform=\"rotate(0, 314.269, 397.6)\" x=\"314.269\" y=\"397.6\">p</text>\n",
"</g>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:16; text-anchor:middle;\" transform=\"rotate(-90, 14.4, 195.808)\" x=\"14.4\" y=\"195.808\">quantity</text>\n",
"</g>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#009af9; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,31.4961 49.5416,36.6206 50.8607,41.5858 52.1797,46.3999 58.278,66.894 64.3764,84.9189 70.4747,100.934 76.573,115.286 88.7697,140.023 100.966,160.679 \n",
" 115.626,181.51 130.285,199.058 142.679,211.92 155.074,223.332 182.001,244.265 206.826,260.05 233.44,274.25 264.276,288.074 291.679,298.554 315.709,306.644 \n",
" 341.908,314.511 368.263,321.593 397.402,328.614 425.127,334.633 453.954,340.314 478.221,344.704 507.847,349.642 535.697,353.914 580.315,360.121 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#ff0000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,342.737 580.315,342.737 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1002)\" style=\"stroke:#008000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 48.2225,169.619 580.315,169.619 \n",
" \"/>\n",
"<polygon clip-path=\"url(#clip1000)\" points=\"\n",
"445.782,97.7361 562.315,97.7361 562.315,52.3761 445.782,52.3761 \n",
" \" fill=\"#ffffff\" fill-opacity=\"1\"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#000000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 445.782,97.7361 562.315,97.7361 562.315,52.3761 445.782,52.3761 445.782,97.7361 \n",
" \"/>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#009af9; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 451.782,67.4961 487.782,67.4961 \n",
" \"/>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:start;\" transform=\"rotate(0, 493.782, 71.9961)\" x=\"493.782\" y=\"71.9961\">q(p)</text>\n",
"</g>\n",
"<polyline clip-path=\"url(#clip1000)\" style=\"stroke:#ff0000; stroke-width:1; stroke-opacity:1; fill:none\" points=\"\n",
" 451.782,82.6161 487.782,82.6161 \n",
" \"/>\n",
"<g clip-path=\"url(#clip1000)\">\n",
"<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:12; text-anchor:start;\" transform=\"rotate(0, 493.782, 87.1161)\" x=\"493.782\" y=\"87.1161\">qty at p= 4</text>\n",
"</g>\n",
"</svg>\n"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Plot the demand function\n",
"using Plots\n",
"q(p) = 2*p^(-1/2)\n",
"\n",
"plot(q,0.5,5,label=\"q(p)\", title=\"Demand curve\", xaxis=\"p\", yaxis=\"quantity\")\n",
"hline!([1, 2], labels = [\"qty at p= 4\" \"qty at p'= 1\"], color = [\"red\", \"green\"])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Analytical solution"
]
},
{
"cell_type": "code",
"execution_count": 92,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Change in consumer surplus"
]
},
{
"data": {
"text/latex": [
"$ \\Delta CS = \\int_1^{\\infty} q(p)\\, \\mathrm{d}p - \\int_4^{\\infty} q(p)\\, \\mathrm{d}p \n",
"= \\int_1^{4} q(p)\\,\\mathrm{d}p = \\int_1^{4} 2p^{-1/2}\\, \\mathrm{d}p\n",
"= \\int_1^{4} q(p)\\,\\mathrm{d}p = \\int_1^{4} 2p^{-1/2}\\, \\mathrm{d}p\n",
"= [4p^{1/2}]_1^4\n",
"= 4(4^{1/2}-1^{1/2})\n",
"= 4\n",
"$"
],
"text/plain": [
"L\"$ \\Delta CS = \\int_1^{\\infty} q(p)\\, \\mathrm{d}p - \\int_4^{\\infty} q(p)\\, \\mathrm{d}p \n",
"= \\int_1^{4} q(p)\\,\\mathrm{d}p = \\int_1^{4} 2p^{-1/2}\\, \\mathrm{d}p\n",
"= \\int_1^{4} q(p)\\,\\mathrm{d}p = \\int_1^{4} 2p^{-1/2}\\, \\mathrm{d}p\n",
"= [4p^{1/2}]_1^4\n",
"= 4(4^{1/2}-1^{1/2})\n",
"= 4\n",
"$\""
]
},
"execution_count": 92,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using LaTeXStrings\n",
"\n",
"print(\"Change in consumer surplus\")\n",
"L\"\"\" \\Delta CS = \\int_1^{\\infty} q(p)\\, \\mathrm{d}p - \\int_4^{\\infty} q(p)\\, \\mathrm{d}p \n",
"= \\int_1^{4} q(p)\\,\\mathrm{d}p = \\int_1^{4} 2p^{-1/2}\\, \\mathrm{d}p\n",
"= \\int_1^{4} q(p)\\,\\mathrm{d}p = \\int_1^{4} 2p^{-1/2}\\, \\mathrm{d}p\n",
"= [4p^{1/2}]_1^4\n",
"= 4(4^{1/2}-1^{1/2})\n",
"= 4\n",
"\"\"\"\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 0.6.2",
"language": "julia",
"name": "julia-0.6"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "0.6.2"
}
},
"nbformat": 4,
"nbformat_minor": 1
}
Loading