{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Homework 4 Solutions\n", "\n", "## Problem 1 (30 points)\n", "\n", "### Problem 1.1 (8 points)\n", "\n", "With $Y = B_0(1-R)$, $B_0 = 50$ a constant, $\\mathbb{E}[R] = 0.6$ and\n", "$\\text{Var}(R) = 0.01$:\n", "\n", "$$\\mathbb{E}[Y] = B_0\\left(1 - \\mathbb{E}[R]\\right) = 50(1 - 0.6) = 20\\ \\text{mg/L}$$\n", "\n", "$$\\text{Var}(Y) = B_0^2\\,\\text{Var}(R) = 50^2 (0.01) = 25\\ (\\text{mg/L})^2, \\qquad \\text{so } \\sigma_Y = 5\\ \\text{mg/L}.$$\n", "\n", "For $Z = Y^2$, use the identity\n", "$\\text{Var}(Y) = \\mathbb{E}[Y^2] - \\left(\\mathbb{E}[Y]\\right)^2$,\n", "rearranged:\n", "\n", "$$\\mathbb{E}[Z] = \\mathbb{E}[Y^2] = \\text{Var}(Y) + \\left(\\mathbb{E}[Y]\\right)^2 = 25 + 400 = 425\\ (\\text{mg/L})^2,$$\n", "\n", "against $\\left(\\mathbb{E}[Y]\\right)^2 = 400\\ (\\text{mg/L})^2$. They\n", "differ by exactly $25\\ (\\text{mg/L})^2$.\n", "\n", "This implies that running a model once at the mean input gives\n", "$h(\\mathbb{E}[X])$, but the quantity you almost always want is the\n", "average outcome $\\mathbb{E}[h(X)]$. These are equal only when $h$ is\n", "linear. For any curved model the single run at the mean is biased.\n", "\n", "### Problem 1.2 (7 points)\n", "\n", "The 95% confidence interval is $0.08 \\pm 1.96(0.006) = 0.08 \\pm 0.0118$,\n", "so $[0.068, 0.092]$.\n", "\n", "The standard error scales as $\\sigma_Y/\\sqrt{n}$, so reducing it by a\n", "factor $f$ costs $f^2$ times as many samples. Here\n", "$f = 0.006/0.002 = 3$, so\n", "\n", "$$n_\\text{new} = 3^2 \\times 2\\text{,}000 = 18\\text{,}000\\ \\text{samples}.$$\n", "\n", "> **Note**\n", ">\n", "> Estimating $\\sigma_Y$ from $\\hat{p}$ instead, as\n", "> $\\sqrt{\\hat{p}(1-\\hat{p})}$, gives\n", "> $n = 0.08(0.92)/0.002^2 = 18\\text{,}400$, which agrees with\n", "> $18\\text{,}000$ to the precision of the rounded standard error\n", "> $0.006$.\n", "\n", "A finer spatial grid will not help, as it reduces *discretization*\n", "error. The $0.006$ here is *Monte Carlo sampling* error, which comes\n", "from having drawn a finite number of random inputs.\n", "\n", "### Problem 1.3 (7 points)\n", "\n", "The annual exceedance probability corresponding to a 50-year return\n", "period is $p = 1/50 = 0.02$, or 2% per year.\n", "\n", "Assuming independent years, the probability of *no* exceedance is\n", "$(1-p)^{20}$, so\n", "\n", "$$\\mathbb{P}(\\text{at least one}) = 1 - (1 - 0.02)^{20} \\approx 0.33.$$\n", "\n", "The problem with the client’s comment is that a return period is a\n", "long-run **average** frequency. Last year’s flood carries no information\n", "about this year’s due to the assumption about independence. The chance\n", "of overtopping next year, or any subsequent year, is still 2%.\n", "\n", "### Problem 1.4 (8 points)\n", "\n", "| | Uncertainty | Type |\n", "|:----|:---------------------------------------------|:--------------|\n", "| 1 | Daily variation in influent flow | **Aleatory** |\n", "| 2 | Strength of the ice-albedo feedback | **Epistemic** |\n", "| 3 | Rate of sediment recycling in a shallow lake | **Epistemic** |\n", "| 4 | Which of next year’s storms is largest | **Aleatory** |\n", "\n", "1. The flow genuinely varies day to day. Collecting more data could\n", " help us characterise that variability better but does not\n", " necessarily reduce or remove it.\n", "2. The feedback has a definite strength, set by the physics of the\n", " climate system; what is uncertain is our knowledge of it.\n", "3. For a given lake, recycling proceeds at a rate set by its sediments\n", " and chemistry. We are uncertain because it is hard to measure, not\n", " because of true randomness in the rate.\n", "4. Which storm will turn out to be the largest is intrinsically random.\n", "\n", "## Problem 2 (20 points)\n", "\n", "### Problem 2.1 (5 points)" ], "id": "ee5e6957-de32-4ed0-a3bf-7e9fd3b7d6fa" }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "loss_limit, steepness, midpoint = 200_000, 0.8, 3\n", "\n", "function depth_damage(depth)\n", " if depth <= 0\n", " return 0\n", " end\n", " return loss_limit / (1 + exp(-steepness * (depth - midpoint)))\n", "end\n", "\n", "hazard = LogNormal(1.2, 0.3)" ], "id": "4" }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [], "id": "cell-fig-dd" }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# Change in damage per metre of depth, measured over a 2 cm step around the midpoint\n", "peak_slope = (depth_damage(midpoint + 0.01) - depth_damage(midpoint - 0.01)) / 0.02" ], "id": "8" }, { "cell_type": "markdown", "metadata": {}, "source": [ "The damage curve is steepest around its midpoint, $h_0 = 3$ m, where the\n", "logistic curve is most nearly vertical — roughly the 2 m to 4 m range.\n", "At the midpoint the damage rises by about \\$40000 per metre of depth.\n", "That band is also where most of the hazard distribution’s mass sits, so\n", "the structure is squarely in the region where a small change in depth\n", "produces a large change in loss. Away from that band, in shallow water\n", "or as the house nears a total loss, additional depth still adds damage,\n", "but more slowly.\n", "\n", "### Problem 2.2 (8 points)" ], "id": "0d994b7c-adeb-4f8a-9f16-95e5ee3ffcc2" }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "Random.seed!(4750)\n", "\n", "function damage_sample(n)\n", " losses = zeros(n)\n", " for i in 1:n\n", " losses[i] = depth_damage(rand(hazard))\n", " end\n", " return losses\n", "end\n", "\n", "sample_sizes = [1_000, 2_000, 4_000, 8_000, 16_000, 32_000]\n", "estimates = zeros(length(sample_sizes))\n", "half_widths = zeros(length(sample_sizes))\n", "for i in 1:length(sample_sizes)\n", " losses = damage_sample(sample_sizes[i])\n", " estimates[i] = mean(losses)\n", " half_widths[i] = 1.96 * std(losses) / sqrt(sample_sizes[i])\n", "end" ], "id": "10" }, { "cell_type": "markdown", "metadata": {}, "source": [ "| $n$ | estimate | 95% half-width | change on doubling |\n", "|--------:|---------:|---------------:|-------------------:|\n", "| $1000$ | \\$114724 | ±\\$2115 | — |\n", "| $2000$ | \\$114684 | ±\\$1515 | \\$39 |\n", "| $4000$ | \\$114931 | ±\\$1063 | \\$247 |\n", "| $8000$ | \\$114658 | ±\\$743 | \\$273 |\n", "| $16000$ | \\$115279 | ±\\$529 | \\$621 |\n", "| $32000$ | \\$115541 | ±\\$378 | \\$263 |\n", "\n", "At $n = 32\\text{,}000$ the expected annual damage is \\$115541, with a\n", "95% confidence interval of roughly ±\\$378. Seed `4750`.\n", "\n", "The quantity being estimated is the mean annual damage. Across the whole\n", "table the estimate never moves by more than a few hundred dollars on a\n", "six-figure loss — under 1% — and the reported half-width shrinks\n", "steadily as $1/\\sqrt{n}$. By $n = 8\\text{,}000$ the estimate is stable\n", "well inside any precision this decision needs, and the extra doublings\n", "past that change nothing of consequence.\n", "\n", "### Problem 2.3 (7 points)" ], "id": "caaedff1-9913-4e76-bec3-0fca83495742" }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "Random.seed!(20)\n", "big_sample = damage_sample(200_000)\n", "mean_damage = mean(big_sample)\n", "tail_damage = quantile(big_sample, 0.99)" ], "id": "14" }, { "cell_type": "markdown", "metadata": {}, "source": [ "The mean annual damage is \\$115086, while the 99th percentile is\n", "\\$189819 — about \\$74733 more, or 1.6 times the mean, and close to the\n", "\\$200,000 ceiling of the damage function.\n", "\n", "Why is a percentile less precise? The estimate of the mean uses each of\n", "the $n$ draws. In contrast, the 99th percentile is determined by the\n", "draws that land above or near it. Only 1% of the sample lies above it,\n", "so draws in that part of the distribution are sparse, and exactly where\n", "those few happen to fall moves the estimate.\n", "\n", "Flood depth only captures the hazard, it cannot say where the properties\n", "are that would be exposed to flooding or how vulnerable those properties\n", "are to a certain level of flooding. The depth-damage function is the\n", "**vulnerability** component: it converts a physical hazard (a depth)\n", "into a consequence (a loss)." ], "id": "4568cfeb-b275-4eb4-a8ba-0589cd313937" } ], "nbformat": 4, "nbformat_minor": 5, "metadata": { "kernel_info": { "name": "julia" }, "kernelspec": { "name": "julia", "display_name": "Julia", "language": "julia" }, "language_info": { "name": "julia", "codemirror_mode": "julia", "version": "1.11.5" } } }