From 4ea4668f3cffe95e48549464e2718d8d624860e4 Mon Sep 17 00:00:00 2001 From: John Stachurski Date: Fri, 3 Mar 2023 17:11:23 +1100 Subject: [PATCH 1/9] misc --- cross_section.md | 647 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 647 insertions(+) create mode 100644 cross_section.md diff --git a/cross_section.md b/cross_section.md new file mode 100644 index 000000000..3c061810c --- /dev/null +++ b/cross_section.md @@ -0,0 +1,647 @@ +--- +jupytext: + text_representation: + extension: .md + format_name: myst +kernelspec: + display_name: Python 3 + language: python + name: python3 +--- + + +# Heavy-Tailed Distributions + +```{contents} Contents +:depth: 2 +``` + +In addition to what's in Anaconda, this lecture will need the following libraries: + +```{code-cell} ipython +--- +tags: [hide-output] +--- +!pip install quantecon +!pip install --upgrade yfinance +``` +We run the following code to prepare for the lecture: + +```{code-cell} ipython +%matplotlib inline +import matplotlib.pyplot as plt +plt.rcParams["figure.figsize"] = (11, 5) #set default figure size +import numpy as np +import quantecon as qe +from scipy.stats import norm +import yfinance as yf +import pandas as pd +from pandas.plotting import register_matplotlib_converters +register_matplotlib_converters() +``` + +## Overview + +TODO -- explain that there is a more detailed and more mathematical discussion of heavy tails in https://python.quantecon.org/heavy_tails.html + +In this section we give some motivation for the lecture. + +### Introduction: Light Tails + +Most commonly used probability distributions in classical statistics and +the natural sciences have "light tails." + +To explain this concept, let's look first at examples. + +The classic example is the [normal distribution](https://en.wikipedia.org/wiki/Normal_distribution), which has density + +$$ f(x) = \frac{1}{\sqrt{2\pi}\sigma} + \exp\left( -\frac{(x-\mu)^2}{2 \sigma^2} \right) +$$ + +on the real line $\mathbb R = (-\infty, \infty)$. + +The two parameters $\mu$ and $\sigma$ are the mean and standard deviation +respectively. + +As $x$ deviates from $\mu$, the value of $f(x)$ goes to zero extremely +quickly. + +We can see this when we plot the density and show a histogram of observations, +as with the following code (which assumes $\mu=0$ and $\sigma=1$). + +```{code-cell} ipython +fig, ax = plt.subplots() +X = norm.rvs(size=1_000_000) +ax.hist(X, bins=40, alpha=0.4, label='histogram', density=True) +x_grid = np.linspace(-4, 4, 400) +ax.plot(x_grid, norm.pdf(x_grid), label='density') +ax.legend() +plt.show() +``` + +Notice how + +* the density's tails converge quickly to zero in both directions and +* even with 1,000,000 draws, we get no very large or very small observations. + +We can see the last point more clearly by executing + +```{code-cell} ipython +X.min(), X.max() +``` + +Here's another view of draws from the same distribution: + +```{code-cell} python3 +n = 2000 +fig, ax = plt.subplots() +data = norm.rvs(size=n) +ax.plot(list(range(n)), data, linestyle='', marker='o', alpha=0.5, ms=4) +ax.vlines(list(range(n)), 0, data, lw=0.2) +ax.set_ylim(-15, 15) +ax.set_xlabel('$i$') +ax.set_ylabel('$X_i$', rotation=0) +plt.show() +``` + +We have plotted each individual draw $X_i$ against $i$. + +None are very large or very small. + +In other words, extreme observations are rare and draws tend not to deviate +too much from the mean. + +As a result, many statisticians and econometricians +use rules of thumb such as "outcomes more than four or five +standard deviations from the mean can safely be ignored." + + +### When Are Light Tails Valid? + +Distributions that rarely generate extreme values are called light-tailed. + +For example, human height is light-tailed. + +Yes, it's true that we see some very tall people. + +* For example, basketballer [Sun Mingming](https://en.wikipedia.org/wiki/Sun_Mingming) is 2.32 meters tall + +But have you ever heard of someone who is 20 meters tall? Or 200? Or 2000? + +Have you ever wondered why not? + +After all, there are 8 billion people in the world! + +In essence, the reason we don't see such draws is that the distribution of +human high has very light tails. + +In fact human height is approximately normally distributed. + + +### Returns on Assets + + +But now we have to ask: does economic data always look like this? + +Let's look at some financial data first. + +Our aim is to plot the daily change in the price of Amazon (AMZN) stock for +the period from 1st January 2015 to 1st July 2022. + +This equates to daily returns if we set dividends aside. + +The code below produces the desired plot using Yahoo financial data via the `yfinance` library. + +```{code-cell} python3 +s = yf.download('AMZN', '2015-1-1', '2022-7-1')['Adj Close'] +r = s.pct_change() + +fig, ax = plt.subplots() + +ax.plot(r, linestyle='', marker='o', alpha=0.5, ms=4) +ax.vlines(r.index, 0, r.values, lw=0.2) +ax.set_ylabel('returns', fontsize=12) +ax.set_xlabel('date', fontsize=12) + +plt.show() +``` + +This data looks different to the draws from the normal distribution. + +Several of observations are quite extreme. + +We get a similar picture if we look at other assets, such as Bitcoin + +```{code-cell} python3 +s = yf.download('BTC-USD', '2015-1-1', '2022-7-1')['Adj Close'] +r = s.pct_change() + +fig, ax = plt.subplots() + +ax.plot(r, linestyle='', marker='o', alpha=0.5, ms=4) +ax.vlines(r.index, 0, r.values, lw=0.2) +ax.set_ylabel('returns', fontsize=12) +ax.set_xlabel('date', fontsize=12) + +plt.show() +``` + +The histogram also looks different to the histogram of the normal +distribution: + + +```{code-cell} python3 +fig, ax = plt.subplots() +ax.hist(r, bins=60, alpha=0.4, label='bitcoin returns', density=True) +ax.set_xlabel('returns', fontsize=12) +plt.show() +``` + +If we look at higher frequency returns data (e.g., tick-by-tick), we often see even more +extreme observations. + +See, for example, {cite}`mandelbrot1963variation` or {cite}`rachev2003handbook`. + + +### Other Data + +The data we have just seen is said to be "heavy-tailed". + +With heavy-tailed distributions, extreme outcomes occur relatively +frequently. + +(A more careful definition is given below) + +Importantly, there are many examples of heavy-tailed distributions +observed in economic and financial settings include + +For example, the income and the wealth distributions are heavy-tailed (see, e.g., {cite}`pareto1896cours`, {cite}`benhabib2018skewed`). + +* You can imagine this: most people have low or modest wealth but some people + are extremely rich. + +The firm size distribution is also heavy-tailed ({cite}`axtell2001zipf`, {cite}`gabaix2016power`}). + +* You can imagine this too: most firms are small but some firms are enormous. + +The distribution of town and city sizes is heavy-tailed ({cite}`rozenfeld2011area`, {cite}`gabaix2016power`). + +* Most towns and cities are small but some are very large. + + +### Why Should We Care? + +Heavy tails are common in economic data but does that mean they are important? + +The answer to this question is affirmative! + +When distributions are heavy-tailed, we need to think carefully about issues +like + +* diversification and risk +* forecasting +* taxation (across a heavy-tailed income distribution), etc. + +We return to these points below. + + + +## Visual Comparisons + +Let's do some more visual comparisons to help us build intuition on the +difference between light and heavy tails. + + +The figure below shows a simulation. (You will be asked to replicate it in +the exercises.) + +The top two subfigures each show 120 independent draws from the normal +distribution, which is light-tailed. + +The bottom subfigure shows 120 independent draws from [the Cauchy +distribution](https://en.wikipedia.org/wiki/Cauchy_distribution), which is +heavy-tailed. + +(light_heavy_fig1)= +```{figure} /_static/lecture_specific/heavy_tails/light_heavy_fig1.png + +``` + +In the top subfigure, the standard deviation of the normal distribution is 2, +and the draws are clustered around the mean. + +In the middle subfigure, the standard deviation is increased to 12 and, as +expected, the amount of dispersion rises. + +The bottom subfigure, with the Cauchy draws, shows a different pattern: tight +clustering around the mean for the great majority of observations, combined +with a few sudden large deviations from the mean. + +This is typical of a heavy-tailed distribution. + + +## Heavy Tails in Economic Cross-Sectional Distributions + +TODO + +- Shu please add data and plots. +- Please use empirical CCDF, as in the macro dynamics book + + + +TODO Review exercises below --- are they all too hard for undergrads or should +we keep some of them. + + +## Exercises + +```{exercise} +:label: ht_ex1 + +Replicate {ref}`the figure presented above ` that compares normal and Cauchy draws. + +Use `np.random.seed(11)` to set the seed. +``` + + +```{exercise} +:label: ht_ex2 + +Prove: If $X$ has a Pareto tail with tail index $\alpha$, then +$\mathbb E[X^r] = \infty$ for all $r \geq \alpha$. +``` + + +```{exercise} +:label: ht_ex3 + +Repeat exercise 1, but replace the three distributions (two normal, one +Cauchy) with three Pareto distributions using different choices of +$\alpha$. + +For $\alpha$, try 1.15, 1.5 and 1.75. + +Use `np.random.seed(11)` to set the seed. +``` + + +```{exercise} +:label: ht_ex4 + +Replicate the rank-size plot figure {ref}`presented above `. + +If you like you can use the function `qe.rank_size` from the `quantecon` library to generate the plots. + +Use `np.random.seed(13)` to set the seed. +``` + +```{exercise} +:label: ht_ex5 + +There is an ongoing argument about whether the firm size distribution should +be modeled as a Pareto distribution or a lognormal distribution (see, e.g., +{cite}`fujiwara2004pareto`, {cite}`kondo2018us` or {cite}`schluter2019size`). + +This sounds esoteric but has real implications for a variety of economic +phenomena. + +To illustrate this fact in a simple way, let us consider an economy with +100,000 firms, an interest rate of `r = 0.05` and a corporate tax rate of +15%. + +Your task is to estimate the present discounted value of projected corporate +tax revenue over the next 10 years. + +Because we are forecasting, we need a model. + +We will suppose that + +1. the number of firms and the firm size distribution (measured in profits) remain fixed and +1. the firm size distribution is either lognormal or Pareto. + +Present discounted value of tax revenue will be estimated by + +1. generating 100,000 draws of firm profit from the firm size distribution, +1. multiplying by the tax rate, and +1. summing the results with discounting to obtain present value. + +The Pareto distribution is assumed to take the form {eq}`pareto` with $\bar x = 1$ and $\alpha = 1.05$. + +(The value the tail index $\alpha$ is plausible given the data {cite}`gabaix2016power`.) + +To make the lognormal option as similar as possible to the Pareto option, choose its parameters such that the mean and median of both distributions are the same. + +Note that, for each distribution, your estimate of tax revenue will be random because it is based on a finite number of draws. + +To take this into account, generate 100 replications (evaluations of tax revenue) for each of the two distributions and compare the two samples by + +* producing a [violin plot](https://en.wikipedia.org/wiki/Violin_plot) visualizing the two samples side-by-side and +* printing the mean and standard deviation of both samples. + +For the seed use `np.random.seed(1234)`. + +What differences do you observe? + +(Note: a better approach to this problem would be to model firm dynamics and +try to track individual firms given the current distribution. We will discuss +firm dynamics in later lectures.) +``` + +## Solutions + +```{solution-start} ht_ex1 +:class: dropdown +``` + +```{code-cell} python3 +n = 120 +np.random.seed(11) + +fig, axes = plt.subplots(3, 1, figsize=(6, 12)) + +for ax in axes: + ax.set_ylim((-120, 120)) + +s_vals = 2, 12 + +for ax, s in zip(axes[:2], s_vals): + data = np.random.randn(n) * s + ax.plot(list(range(n)), data, linestyle='', marker='o', alpha=0.5, ms=4) + ax.vlines(list(range(n)), 0, data, lw=0.2) + ax.set_title(f"draws from $N(0, \sigma^2)$ with $\sigma = {s}$", fontsize=11) + +ax = axes[2] +distribution = cauchy() +data = distribution.rvs(n) +ax.plot(list(range(n)), data, linestyle='', marker='o', alpha=0.5, ms=4) +ax.vlines(list(range(n)), 0, data, lw=0.2) +ax.set_title(f"draws from the Cauchy distribution", fontsize=11) + +plt.subplots_adjust(hspace=0.25) + +plt.show() +``` + +```{solution-end} +``` + + +```{solution-start} ht_ex2 +:class: dropdown + +Let $X$ have a Pareto tail with tail index $\alpha$ and let $F$ be its cdf. + +Fix $r \geq \alpha$. + +As discussed after {eq}`plrt`, we can take positive constants $b$ and $\bar x$ such that + +$$ +\mathbb P\{X > x\} \geq b x^{- \alpha} \text{ whenever } x \geq \bar x +$$ + +But then + +$$ +\mathbb E X^r = r \int_0^\infty x^{r-1} \mathbb P\{ X > x \} x +\geq +r \int_0^{\bar x} x^{r-1} \mathbb P\{ X > x \} x ++ r \int_{\bar x}^\infty x^{r-1} b x^{-\alpha} x. +$$ + +We know that $\int_{\bar x}^\infty x^{r-\alpha-1} x = \infty$ whenever $r - \alpha - 1 \geq -1$. + +Since $r \geq \alpha$, we have $\mathbb E X^r = \infty$. +``` + + +```{solution-end} +``` + +```{solution-start} ht_ex3 +:class: dropdown +``` + +```{code-cell} ipython3 +from scipy.stats import pareto + +np.random.seed(11) + +n = 120 +alphas = [1.15, 1.50, 1.75] + +fig, axes = plt.subplots(3, 1, figsize=(6, 8)) + +for (a, ax) in zip(alphas, axes): + ax.set_ylim((-5, 50)) + data = pareto.rvs(size=n, scale=1, b=a) + ax.plot(list(range(n)), data, linestyle='', marker='o', alpha=0.5, ms=4) + ax.vlines(list(range(n)), 0, data, lw=0.2) + ax.set_title(f"Pareto draws with $\\alpha = {a}$", fontsize=11) + +plt.subplots_adjust(hspace=0.4) + +plt.show() +``` + +```{solution-end} +``` + + +```{solution-start} ht_ex4 +:class: dropdown +``` + +First let's generate the data for the plots: + +```{code-cell} ipython3 +sample_size = 1000 +np.random.seed(13) +z = np.random.randn(sample_size) + +data_1 = np.abs(z) +data_2 = np.exp(z) +data_3 = np.exp(np.random.exponential(scale=1.0, size=sample_size)) + +data_list = [data_1, data_2, data_3] +``` + +Now we plot the data: + +```{code-cell} ipython3 +fig, axes = plt.subplots(3, 1, figsize=(6, 8)) +axes = axes.flatten() +labels = ['$|z|$', '$\exp(z)$', 'Pareto with tail index $1.0$'] + +for data, label, ax in zip(data_list, labels, axes): + + rank_data, size_data = qe.rank_size(data) + + ax.loglog(rank_data, size_data, 'o', markersize=3.0, alpha=0.5, label=label) + ax.set_xlabel("log rank") + ax.set_ylabel("log size") + + ax.legend() + +fig.subplots_adjust(hspace=0.4) + +plt.show() +``` + +```{solution-end} +``` + + + +```{solution-start} ht_ex5 +:class: dropdown +``` + +To do the exercise, we need to choose the parameters $\mu$ +and $\sigma$ of the lognormal distribution to match the mean and median +of the Pareto distribution. + +Here we understand the lognormal distribution as that of the random variable +$\exp(\mu + \sigma Z)$ when $Z$ is standard normal. + +The mean and median of the Pareto distribution {eq}`pareto` with +$\bar x = 1$ are + +$$ +\text{mean } = \frac{\alpha}{\alpha - 1} +\quad \text{and} \quad +\text{median } = 2^{1/\alpha} +$$ + +Using the corresponding expressions for the lognormal distribution leads us to +the equations + +$$ +\frac{\alpha}{\alpha - 1} = \exp(\mu + \sigma^2/2) +\quad \text{and} \quad +2^{1/\alpha} = \exp(\mu) +$$ + +which we solve for $\mu$ and $\sigma$ given $\alpha = 1.05$. + +Here is code that generates the two samples, produces the violin plot and +prints the mean and standard deviation of the two samples. + +```{code-cell} ipython3 +num_firms = 100_000 +num_years = 10 +tax_rate = 0.15 +r = 0.05 + +β = 1 / (1 + r) # discount factor + +x_bar = 1.0 +α = 1.05 + +def pareto_rvs(n): + "Uses a standard method to generate Pareto draws." + u = np.random.uniform(size=n) + y = x_bar / (u**(1/α)) + return y +``` + +Let's compute the lognormal parameters: + +```{code-cell} ipython3 +μ = np.log(2) / α +σ_sq = 2 * (np.log(α/(α - 1)) - np.log(2)/α) +σ = np.sqrt(σ_sq) +``` + +Here's a function to compute a single estimate of tax revenue for a particular +choice of distribution `dist`. + +```{code-cell} ipython3 +def tax_rev(dist): + tax_raised = 0 + for t in range(num_years): + if dist == 'pareto': + π = pareto_rvs(num_firms) + else: + π = np.exp(μ + σ * np.random.randn(num_firms)) + tax_raised += β**t * np.sum(π * tax_rate) + return tax_raised +``` + +Now let's generate the violin plot. + +```{code-cell} ipython3 +num_reps = 100 +np.random.seed(1234) + +tax_rev_lognorm = np.empty(num_reps) +tax_rev_pareto = np.empty(num_reps) + +for i in range(num_reps): + tax_rev_pareto[i] = tax_rev('pareto') + tax_rev_lognorm[i] = tax_rev('lognorm') + +fig, ax = plt.subplots() + +data = tax_rev_pareto, tax_rev_lognorm + +ax.violinplot(data) + +plt.show() +``` + +Finally, let's print the means and standard deviations. + +```{code-cell} ipython3 +tax_rev_pareto.mean(), tax_rev_pareto.std() +``` + +```{code-cell} ipython3 +tax_rev_lognorm.mean(), tax_rev_lognorm.std() +``` + +Looking at the output of the code, our main conclusion is that the Pareto +assumption leads to a lower mean and greater dispersion. + +```{solution-end} +``` From 131da2172a8793b209188d2c3e5b364b706c0264 Mon Sep 17 00:00:00 2001 From: Shu Date: Sun, 16 Apr 2023 16:27:12 +1000 Subject: [PATCH 2/9] embed_figure_code_w_fix --- cross_section.md | 199 +++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 183 insertions(+), 16 deletions(-) diff --git a/cross_section.md b/cross_section.md index 3c061810c..11754b665 100644 --- a/cross_section.md +++ b/cross_section.md @@ -3,13 +3,14 @@ jupytext: text_representation: extension: .md format_name: myst + format_version: 0.13 + jupytext_version: 1.14.5 kernelspec: - display_name: Python 3 + display_name: Python 3 (ipykernel) language: python name: python3 --- - # Heavy-Tailed Distributions ```{contents} Contents @@ -18,24 +19,30 @@ kernelspec: In addition to what's in Anaconda, this lecture will need the following libraries: -```{code-cell} ipython ---- -tags: [hide-output] ---- +```{code-cell} ipython3 +:tags: [hide-output] + !pip install quantecon !pip install --upgrade yfinance +!pip install pandas_datareader ``` + We run the following code to prepare for the lecture: -```{code-cell} ipython +```{code-cell} ipython3 %matplotlib inline import matplotlib.pyplot as plt plt.rcParams["figure.figsize"] = (11, 5) #set default figure size import numpy as np import quantecon as qe -from scipy.stats import norm import yfinance as yf import pandas as pd +import pandas_datareader.data as web +import statsmodels.api as sm + +from interpolation import interp +from pandas_datareader import wb +from scipy.stats import norm, cauchy from pandas.plotting import register_matplotlib_converters register_matplotlib_converters() ``` @@ -70,7 +77,7 @@ quickly. We can see this when we plot the density and show a histogram of observations, as with the following code (which assumes $\mu=0$ and $\sigma=1$). -```{code-cell} ipython +```{code-cell} ipython3 fig, ax = plt.subplots() X = norm.rvs(size=1_000_000) ax.hist(X, bins=40, alpha=0.4, label='histogram', density=True) @@ -87,13 +94,13 @@ Notice how We can see the last point more clearly by executing -```{code-cell} ipython +```{code-cell} ipython3 X.min(), X.max() ``` Here's another view of draws from the same distribution: -```{code-cell} python3 +```{code-cell} ipython3 n = 2000 fig, ax = plt.subplots() data = norm.rvs(size=n) @@ -153,7 +160,7 @@ This equates to daily returns if we set dividends aside. The code below produces the desired plot using Yahoo financial data via the `yfinance` library. -```{code-cell} python3 +```{code-cell} ipython3 s = yf.download('AMZN', '2015-1-1', '2022-7-1')['Adj Close'] r = s.pct_change() @@ -173,7 +180,7 @@ Several of observations are quite extreme. We get a similar picture if we look at other assets, such as Bitcoin -```{code-cell} python3 +```{code-cell} ipython3 s = yf.download('BTC-USD', '2015-1-1', '2022-7-1')['Adj Close'] r = s.pct_change() @@ -190,8 +197,7 @@ plt.show() The histogram also looks different to the histogram of the normal distribution: - -```{code-cell} python3 +```{code-cell} ipython3 fig, ax = plt.subplots() ax.hist(r, bins=60, alpha=0.4, label='bitcoin returns', density=True) ax.set_xlabel('returns', fontsize=12) @@ -293,6 +299,167 @@ TODO TODO Review exercises below --- are they all too hard for undergrads or should we keep some of them. +```{code-cell} ipython3 +def extract_wb(varlist=['NY.GDP.MKTP.CD'], c='all', s=1900, e=2021): + df = wb.download(indicator=varlist, country=c, start=s, end=e).stack().unstack(0).reset_index() + df = df.drop(['level_1'], axis=1).set_index(['year']).transpose() + return df +``` + +```{code-cell} ipython3 +def empirical_ccdf(data, + ax, + aw=None, # weights + label=None, + xlabel=None, + add_reg_line=False, + title=None): + """ + Take data vector and return prob values for plotting. + Upgraded empirical_ccdf + """ + y_vals = np.empty_like(data, dtype='float64') + p_vals = np.empty_like(data, dtype='float64') + n = len(data) + if aw is None: + for i, d in enumerate(data): + # record fraction of sample above d + y_vals[i] = np.sum(data >= d) / n + p_vals[i] = np.sum(data == d) / n + else: + fw = np.empty_like(aw, dtype='float64') + for i, a in enumerate(aw): + fw[i] = a / np.sum(aw) + pdf = lambda x: interp(data, fw, x) + data = np.sort(data) + j = 0 + for i, d in enumerate(data): + j += pdf(d) + y_vals[i] = 1- j + + x, y = np.log(data), np.log(y_vals) + + results = sm.OLS(y, sm.add_constant(x)).fit() + b, a = results.params + + kwargs = [('alpha', 0.3)] + if label: + kwargs.append(('label', label)) + kwargs = dict(kwargs) + + ax.scatter(x, y, **kwargs) + if add_reg_line: + ax.plot(x, x * a + b, 'k-', alpha=0.6, label=f"slope = ${a: 1.2f}$") + if not xlabel: + xlabel='log value' + ax.set_xlabel(xlabel, fontsize=12) + ax.set_ylabel("log prob.", fontsize=12) + + if label: + ax.legend(loc='lower left', fontsize=12) + + if title: + ax.set_title(title) + + return np.log(data), y_vals, p_vals +``` + +### GDP + +```{code-cell} ipython3 +df_gdp1 = extract_wb(varlist=['NY.GDP.MKTP.CD']) # gdp for all countries from 1960 to 2022 +df_gdp2 = extract_wb(varlist=['NY.GDP.PCAP.CD']) # gdp per capita for all countries from 1960 to 2022 +``` + +```{code-cell} ipython3 +fig, axes = plt.subplots(1, 2, figsize=(8.8, 3.6)) + +empirical_ccdf(np.asarray(df_gdp1['2021'].dropna()), axes[0], add_reg_line=False, label='GDP') +empirical_ccdf(np.asarray(df_gdp2['2021'].dropna()), axes[1], add_reg_line=False, label='GDP per capita') + +plt.show() +``` + +### Firm size + +```{code-cell} ipython3 +df_fs = pd.read_csv('https://media.githubusercontent.com/media/QuantEcon/high_dim_data/update_csdata/cross_section/forbes-global2000.csv') +df_fs = df_fs[['Country', 'Sales', 'Profits', 'Assets', 'Market Value']] +``` + +```{code-cell} ipython3 +fig, ax = plt.subplots(figsize=(6.4, 3.5)) + +label="firm size (market value)" + +d = df_fs.sort_values('Market Value', ascending=False) + +empirical_ccdf(np.asarray(d['Market Value'])[0:500], ax, label=label, add_reg_line=True) + +plt.show() +``` + +### City size + +```{code-cell} ipython3 +df_cs_us = pd.read_csv('https://raw.githubusercontent.com/QuantEcon/high_dim_data/update_csdata/cross_section/cities_us.txt', delimiter="\t", header=None) +df_cs_us = df_cs_us[[0, 3]] +df_cs_us.columns = 'rank', 'pop' +x = np.asarray(df_cs_us['pop']) +citysize = [] +for i in x: + i = i.replace(",", "") + citysize.append(int(i)) +df_cs_us['pop'] = citysize +``` + +```{code-cell} ipython3 +df_cs_br = pd.read_csv('https://media.githubusercontent.com/media/QuantEcon/high_dim_data/update_csdata/cross_section/cities_brazil.csv', delimiter=",", header=None) +df_cs_br.columns = df_cs_br.iloc[0] +df_cs_br = df_cs_br[1:401] +df_cs_br = df_cs_br.astype({"pop2023": float}) +``` + +```{code-cell} ipython3 +fig, axes = plt.subplots(1, 2, figsize=(8.8, 3.6)) + + +empirical_ccdf(np.asarray(df_cs_us['pop']), axes[0], label="US", add_reg_line=True) +empirical_ccdf(np.asarray(df_cs_br['pop2023']), axes[1], label="Brazil", add_reg_line=True) + +plt.show() +``` + +### Wealth + +```{code-cell} ipython3 +df_w = pd.read_csv('https://media.githubusercontent.com/media/QuantEcon/high_dim_data/update_csdata/cross_section/forbes-billionaires.csv') +df_w = df_w[['country', 'realTimeWorth', 'realTimeRank']].dropna() +df_w = df_w.astype({'realTimeRank': int}) +df_w = df_w.sort_values('realTimeRank', ascending=True).copy() +``` + +```{code-cell} ipython3 +countries = ['United States', 'Japan', 'India', 'Italy'] +N = len(countries) + +fig, axs = plt.subplots(2, 2, figsize=(8, 6)) +axs = axs.flatten() + +for i, c in enumerate(countries): + df_w_c = df_w[df_w['country'] == c].reset_index() + z = np.asarray(df_w_c['realTimeWorth']) + # print('number of the global richest 2000 from '+ c, len(z)) + if len(z) <= 500: # cut-off number: top 500 + z = z[0:500] + + empirical_ccdf(z[0:500], axs[i], label=c, xlabel='log wealth', add_reg_line=True) + +fig.tight_layout() + +plt.show() +``` + ## Exercises @@ -394,7 +561,7 @@ firm dynamics in later lectures.) :class: dropdown ``` -```{code-cell} python3 +```{code-cell} ipython3 n = 120 np.random.seed(11) From d21a8ab01d781458721397b89e1d701a2d050db2 Mon Sep 17 00:00:00 2001 From: Shu Date: Mon, 17 Apr 2023 15:45:02 +1000 Subject: [PATCH 3/9] reorganise_add_toc --- lectures/_toc.yml | 1 + cross_section.md => lectures/cross_section.md | 0 2 files changed, 1 insertion(+) rename cross_section.md => lectures/cross_section.md (100%) diff --git a/lectures/_toc.yml b/lectures/_toc.yml index def7c152e..da76ceda7 100644 --- a/lectures/_toc.yml +++ b/lectures/_toc.yml @@ -8,6 +8,7 @@ parts: - file: long_run_growth - file: business_cycle - file: inequality + - file: cross_section - caption: Tools & Techniques numbered: true chapters: diff --git a/cross_section.md b/lectures/cross_section.md similarity index 100% rename from cross_section.md rename to lectures/cross_section.md From 450250e5036a83ddf41812f123c65fc3c3505ef6 Mon Sep 17 00:00:00 2001 From: mmcky Date: Thu, 20 Apr 2023 07:20:27 +1000 Subject: [PATCH 4/9] fix toc --- lectures/_toc.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lectures/_toc.yml b/lectures/_toc.yml index 266e85bd2..f053246e1 100644 --- a/lectures/_toc.yml +++ b/lectures/_toc.yml @@ -19,11 +19,11 @@ parts: - file: business_cycle - file: inflation_unemployment - file: inequality + - file: distributions - file: cross_section - caption: Useful Tools numbered: true chapters: - - file: distributions - file: geom_series - file: linear_equations - file: eigen_I From 116fa63f4ccaedb4594c7f943b268cd3969dd63b Mon Sep 17 00:00:00 2001 From: Shu Date: Mon, 24 Apr 2023 13:50:37 +1000 Subject: [PATCH 5/9] fix_label --- .../cross_section/light_heavy_fig1.png | Bin 0 -> 31968 bytes .../cross_section/rank_size_fig1.png | Bin 0 -> 25360 bytes lectures/cross_section.md | 63 ++++++++++++++++-- 3 files changed, 56 insertions(+), 7 deletions(-) create mode 100644 lectures/_static/lecture_specific/cross_section/light_heavy_fig1.png create mode 100644 lectures/_static/lecture_specific/cross_section/rank_size_fig1.png diff --git a/lectures/_static/lecture_specific/cross_section/light_heavy_fig1.png b/lectures/_static/lecture_specific/cross_section/light_heavy_fig1.png new file mode 100644 index 0000000000000000000000000000000000000000..2035bd18189d463a33285c28aea6fc798d88bec7 GIT binary patch literal 31968 zcmeFZ1yojD+wZ$*6ltVGFhC`hZb=m+R6;r=q*FSSln|4U6hx3lxm9wl+68c{q94FPJ+y*$H!T{mb7uZ5_2Ca(L^4q||L*48IeZGv5Vb`}FSrSoqevBsrif~*_RDN>WA^|E5GCKT>1X^ z^u$-_=Oxydn8$Bo7^)MCG(&|PXXCtf+&p*3=6=YKO5u}`lHO@*8^yW6!a_@gOASB0 zz>t%Xxo_l4O8WX8wSkPxrOW^Sum6vHCL8uH7{M(J^FAXnF6b zEQY~F)=ZF2aCDQkK<2|<|IC{olWrCdE`8lmgpS< z<7(hXiNnlYdXcB+5WAhFN56mnHkkbp)AIfMEjPFQyndhVdjf`sJIez_mQ)!T8NbHH zh~Sc0S?q;{h37oVkXNr?Kl}Mp390Q?mIz;2vDp@LSzc_5;xZomYI6I|om=ki!sIJA7V}vnpMsShM%yLRn)r+EXh_ZrrZn45T;v+ea6YMHq2W?suP z)-^Nvgx!nfpX}4N?PYI1IX-+H6Z2wg%jx&*?9@n^Bfp9FO4;12`1tTopDtCNo*cSv zd}l@9Qdr0vC*&wsKy+6@!5;3g&YU+nIaw;+l)*q&R#vmXln_a(u4YV4O>ODuz(an| z&7~9;I<@G@+`C7po~8LLEbJUM4$duA3bv5sog1fn%@R`s`NnMQ><=a0CMHHK50_YS zcDxoY>A3@=MtJpG0KseM)7e>5q;}PNSEQ&q3u|O>L&UiC-PO{WNbR0}fe|~bYuB&G z>k)-ZKX19^?Cfj=f1Aq8ttOMgUq9OHAV89{v+p4!0RaL0n{?-Jw@;PDw%-jDK6y(b z6}++G@cPXg7Jm{WcswkutnyZaUuE%`8{mO8mO0KfKo4;`I6C%r2=Mb;L_Xp?K6| zwd%v7;=Oz4nS74&(QEI`aQR$?%Z^pLP5BXvtIXZGdzZ}YXB>o6)%*8ZV9*;Egbfy$ zg-$m|SaR0n<|_Y4l+d>tF20A)?5~$gb#ii&TwJUUP1Uxox}2Pn(s;1Fh`waKys8cc z8yaMstNQ4_I4Z;!78Mch)y!-Ya8+LLCz*PyEKaJbs#<@orcu{r^zv+bY}$T@Ij;rh zuoZvk$cWDR@c4=PrdN3}3`E1m#zv#>Z!Rw^e*SzV7U^Ex=HEPB(5NBn*L8`%t)M{u z;>C-ponc#z@zBsvay2zI76Abb&%MbY^g~mzw$-HT>+4eZ%s&$)xX>mB64EtX&#hz; z78d4h`*hwQ(xGyHHkO~+oIf@oF)`8H+FEU{$bShQyo`ZClJM$B>7a(^I0U#^(7?0$ z`ud%HeOF%63$ajg>adlUm;0liv9VE~+uYoIAl*NRhhj za0ziJ6&EgY>GI`C7(ANgPFIN-g#Fz%X9qpUUOV73C$Y7@+B`Ts#6l{osvtb{{23b? z`>^{<{Pf5XKE%w;O)@w*n4Xi<1Y;#=U_cWa4=;1`UNrYD4qGcLb_)v&1O8YH#I%5> zsHiA7JiM2^doI!E<&PhV+uNxMQNKDnpT*oX8yelo$>H2)RQK6g7$NQ3#xWl^b1FP? z8TS&napMN>=_AHV`8xx~5g#in>3t4|$0ueQ6{yx3sfCzXS@Bqv;#3loJYXo}#g~+o z(bUz|;S&-De1H6n>0}KrcGuC-G2f((a>R8y9Py>Itt5ph)zaRMgFMvGs2?7t8$UhS zU5_0t; zOssf?oZn_7Sf0s?qPn_TI`agfX)z6PhEyT{kf(fx}Z)~?Z2qNycX0x?i(lPsHhG=oGzrwJ{ zN~v)zjyH6AdREp$u7t9(GV~uFSggtkGD=EBFo%NBwcVNFLZE>Vrto6kDG(>T(jPa5oUhko3}*_-6X2AQQL?nQZoDTj z>VHxG0E+`-P%MQX|7}Y2R_Xb8-_^sM>nQSNWc4H$noAKViXOYpP$Fj zw__x{^7f3#)1_K0Vo~meMV1YUbKEQ}Z^%}k{ksVtZHW&tJ3ITiGiT1MaN?4_e#og; z85|gh>q1q}Jk4KAMn>FXeS1>fVfAc+Q81+71yr=)}^Ml~}t14CQZ(2yCK7=hOk;^|1F0^%2sNh?L{69|=>6;C-3gXN*3feM#g;p3CT)u3dy(#b%& zUn3*us`MmgP2f`o}I2%tCSgJO*+Ow>~{GG&p#q&Dh=GL7I|U5bn0S7z+yvqVvZd zLhzSOH4Tl6*RI`lq`|?#VY_3>O{F0IR&Cf_~;O6EwQ0B<`*Q7*apYXA= z&jGnV$sU9^5_EQkS{k^q%;?7^ce50eOMK^!kbe(#^nrG2_&1$@w z!QS59b-Pys18Ixos25qS^!xq$F*Y`~LGw$76}B6sq;Ax*=A(~k@gP;%ql$SE9CxX^ zsM|OYdNw$&?!WHN|7V^5zvkWhr4`8Wh#c}dgI-@cL;KpfZ%m(Y$zlekXj|@<^5{E ze2Im{WwQlY%q=aYw(QTrL9%o!C_Q#poE_oOw2X>!)T{0R7*#{#;Cd9YATv$d< z4ijOycu`effAuE|8R=`hwNej*AYvxDhYu+h78moK{@nlZfv;_{gYx|OFCWdgspT?* zp~o>qbMW9dI7o1CaBax>>vmi7s*3T#=iRFqk)-}b2)Q5k2cJL3gvKU)`h-0>IccnJ z3L7~O#`6XY3=HO$mJd!($2?Y)pkZ|UC9X@SnY@F=X(06vc1LJyB^jhtkH*ePV=Mssd8EI)lAiR$g!n0#h zo>7&CHp{&H-+1Tvy&2u;(ne^7i4AyVk)0=$yA+6XK|+MF4x@U7v#<0UQou5dA@F9S zwK%(<#pCdHpck7ojc)tb_SLvbZ2V_9FFuwUa2xIB-xbYc+y9XDwv*+a@A7yJ6S2e* z_Yc#X-dib5ljvq(U-&6gW&@5NDd}4$e!RFe+KkE_B{Jn<_3RC;);~X9tSp53rL@ug zbeAuW(8(vxZX=<|!t%%TvW#lu%t~>KbF;wu`L$haR(#9JWx7EY%O;wK5Lz3%D>gIo zHiiPjls`7#ep-%R)RhFH4s@2Bo*pCQozq?04|%CceQl^k+(~$@s)$X*@x_HcpSh{8 zZz2C6=*(Nzg;zH?>{`@!X}k#T*E@fh=@q`SQ83LXLHAn%18$OcYK{7iSZ|RN%cWXk zp0`*%8<=(nGECnaZ&}tHPcOO`oj-nG%|N_*wek0ha*V^4G4~?Q=pAw6Rag22uWait zR=gQAFNx#14Gza7jD0bm#hS^b4-50TKJTux9wgp*sp`>|?@}o=T zdN4>g_1rEine;67l|R*DF=^O)g!mMT2>rtE{Dj&0#?AQ|rNXC9Mj`r%?RQ3f=nDIc z84l(Mt+6CGi$V z{OwR8`}xu$B2onFWhXR}x5O}gWyM$phA2opz@S=h;W*Rb1C|6G}VQRJxk5=|b*Cb=|UhnhSc<8%N8hWdjC zevr*>J1Z)lgWP>;WBa#R^K@#Y`a{0hUwxKxw|x`aWxWPRT9~m~9o97Uf4CnW^!)v# z1NwhiG+ffhiP8=$6}M@9mA!)V;(^`X(&cx(#r6IDR0L#F)xgB{mW4eO_|Mol9Ws=%b{GxUf0bRmd0ZlmFolVc2N40je?R;~GLROY3qL*?9 zmM+S&{VH){x#5bLK|AbQPSmiypApWa&~e{V@J?pk*6GP)u$Iebgn^Uq0}^dOX1W!9 zxF!))txl)-x=vUI`=IXESa{DK2|S9-=DzWlYQ%!LBKO2B%BsA)@T?MRY9#zis16RW zw>Wg+a!FrO)kaTs#bzdi_^xie!=tv?JU;Wm(E+`UYZ%pi>}V3BT3%UIT*d3xu1Tq?k|MPR+=5jJ z6NtWOf(rv;0KI=fR3(UfocCiv8#WFm)j7FJr};}6AMOPw3ge9ForNw~r(#V)lm=Yf0Yud&dBmp1;&TLHerO8v1U zj1-8M0_*1h?4r6feRG;PV-LkY4sZO^lBU69z(frRs_?|1TJ zk{cM%smaNMVbDmpyzI=(b=LHc?VtMYn55BMab%j3X3ZF0U0udu(%1KN`=XYeIGe+$ zcxFr5(v0iC-xt&1u5Y4NLi7cny?r)woh^>;!$^tgQ1X9*FAwq@))sRO-5ydy0 z?(UKmZoA(~+cFAKqOm)tMN8T(h52jDgfY)6)GQzq<2;iv1~Z0tx;g)s!*Taf2$*a4 z6|R%jKewwxWE6WQjdU+gNjc%4J$b6-6?>*FH=?3Ta7jZhds;$eZgFkB_=5e}mwHwW ze*$f`wsf9kdF_PTOL|$x`@OE~y#4ZAqd>N)B>v+o{!;}Za-bBFaPyEHY`QBFtd|Il ziHrESNsj7S@EB?ZOMDQi5i8D#h(+RasS;|g>sA#n39jR?T+tPhH@@O?=xuGyA^d?5 zpD0)g{w3w*(fMmRDENhN(K$U8Tq_0nTxz=uhY2w{Bq|7p`Bolvz|fS8>Cs&wME`H<#-(_5rie8>08l48eqD4n}f{ zU3D;+&tm`3(9Y)h`0REBQQh|i%6aYbJkg6A!;?1<8V%*IPas^L<&we1OUg_m7P_OM z)byoJgY$Rf$KZ)GWB7!7J5%T~@*vQLzRdMefxP)uX&G~T%6pXDg|R)>F$^@H^s?T3 zt6_MP&pw^0knKwneUk$}aK{FQTe{VfjH3C4S=k_1yC~;la|`h~oLLSKFUPd@9gkbL zjO7JB2$%OB?C(tdiP+)1Nr4}`@gwA%jQ7)PkY>n4d?7c)s7|-8?(%a#$#LR)oJkkH zAe$@CCV1bThWTWiv!W=Is%=&Y%hxuldyLPZ`r`R;lG-HaSdqiwLuQxH$^|sl_bFsY z<%Mmz{IU1y`iBX49$?xz-i^NJ$8w{T>R?vs_s@h38_RYoOzbb7&dGWae0E3PB!!f4 zS!|w%;0ICZjeFxmsf=b?X!P99)`jR9>1F+3oWl!iVVLYk-(1qnz`_%lR>5(EpFBH7Hl|_Y|lUzu`f? zUnd9+P+mw z-K^(Wp6j*Ak7M_=O!eiVq0C*oG9Bq-07&_`w2hBM?&7oltQkjt*h z%SduwyVzy3!mMY8eT|r!d)(~--AQ6hH%ZcFXV!2c`<+fv&I96R&GCY-IhJjWt@q8W zKkGSb3)FbXjpjyL(`URktP)+seXCyk2y^qG^jUM#P4;GN#Q9Xa&)Q411p8KNp=5ht zP=Gp3(MzAj?Lm<-`fS^_d=oFC54n96+fOZ)H5QQ)>dXkJ$gK{nUd>-w^_?GKkr>aF;VL1Ddp?;b<1vxRO? zk=9Q?u_t|}6ZtmJ9w~kOPVOwtxYd`~S#Gz@wP;{My>yAeDmME`04-9%8lkQtER~r$ zRAUuMpp(+cdts|6GNo5AZ)@JEdQ7@6J9FVjIb&h5eY#xWwnrb~=Z~BdeN=m|e`&d_ z&lN4mjy4Y~SU)!B6D=$r99v(VU!OJiB34`V>l7XMgvjL6I5}bxt9hL>A^%uJZ`X!&A=MUK&# z)011jd_Gj=l7?-si_nhu(w$j!Kbl`EA@+OHJFZC^ZaHUe=qq$ac6Ng7 zCn+a~JLWdWuwBKWU3@`Ei010mtEe_crvRJGgm19z%G+tV7%`FG%vcb!GZUq6)R7=2 z=ykwrVQHEF_&ZsgmIhJ|H=#F*AVA(H);me!E34@f+& zYjt?`Cr57Y-oIzQaf2*E{D2LS#OE;(3Q0EStVvJDfAr{4OKa=%Is@unYmlOLX}e*Q zPQ>6~Ud;q%$KHxd4Z+0Sh|G+PyI#3}FESqz|Meo{p9&jcLIxJNT_GW7-$g|kpB(M4T)P6diXY*2c745>N&e9z8Zg^XLBjpdWkk}#;+lH4 zRsfi4wav|VT-AFliISi^z{tGv&m+@?$}jl#wHvx3Vy5`7L3^P2Un8UZr?TDuz$Lq) z@dm))3F+y1=<4S7@zbY<*?%)=guPBrj#0-C_IRMcG?%x46!NS|w6?3OYYptj>DgG@ z;+T_vHHvmt$0G>o_-n2B(HJw`6z;k^>avc$VRvP8Z=;Q?xqJY=!_c$IJ)pp>qc2WR zM<=W|Q$sqfo0*4)^!xYkg6>;v2s&`@+_{7J!kVF`F}k5Y_4N~Hoj5=R z0y_^Ap+1Le7NRI4^M>XMMb_@2jlOVnX-yf zxM@@46Ov|hD4y;`och)o)W*Avy#ObVUeFHt^5u(Wk;Qos4zQS*m@44wOG`<8%M-WcrY2!z zTp7YF`|dhKY2ay1T+(<`cX8NZeM**>mpAPG3mwZ5+=1IITiOOYVQJPFP5nD|+M?to zvW)anFgO08FOMlp9B*~6aQd^BpByZX>S9OHC=X+X_G&STGm{7LMm`#W>g~4<8zUoFC0&j0@%~qmzqEMVbH0CpI=spg@N|Dlj#K5U>Z9 z^5tdk{dQY#KN$a?Qd5~FB$#w{b?+}n{cHKVZYO{QlUQtv1mfQ7goKw}UH4tiYrpL) zKH5r?Xar~C75wu2b@imOvPWLFAmQFtQj)uS_bd^;;KWOoipUFUwzj#q!IYfqNTj;* z_RgoYv{3LmD=IXJGPASKH%Bm0N#S=64P7(+5zP!v1BdrPcVY;HZty1uzM6!(9<4Wz zK2l6OLj6JzlBa%(tWtm@>&J8?xinU#4bCO%Se~A&o|YH55>|N%MUCfHr~K63`ZNn4pISWsKs6X~iF-d?sxKZK z95g^m#l_=u;#u;}jazHjVWJkJv%C8u*q#s|HA-f-{&cg5iRn5yJHH^H;tU3FyK{JW zG18&pEY^=!Sc3`8oj-q`pX#`7NA}jOXM2nKrx`+tY4bH4eR@8}&iSV8w9E~)i-VvB zUlE5v$#Zque3B&j{f7^<4NXlMxyD9DDi0n^H<1yJH-=K#Pc@uj*UUq}1P5mU8v)HY zxw?9+J=O7tviAo0(XP*_KxJhP{pTnyeYxAW@jX2~k(L4P_w%zq-tJc@)6XI-=6YOn?US8$H8?(e_@DxD71p?(zCN&wu~w@OSA6#$4{*JFwH=M+3l0^LxVU&J zAP^WWRUCeXzK{M2fd3NoY}j9Ky1TpEfkizDvUOwY`|B78_#^>(xp)XTwIDB-76cOs zdG2wKd2Zgj1mgboLN9O#C(ja(ezSs~EC)?P+;-4B0D3IXpbitMMJFK4psTxW;7PxV zjt++BwRdn}1peLLR=0BM2v;BOFx6K*@~0~!FCgyG`Rupp_q=-JtQ`o?(IwC!>3FVMkr_-h{MCdY zBQw(;LfFz^0UivvCjFBQsuB0a>^<op^%Hj4gh+bQ3k1rNa{}V6}8T*`Hg7Ws&8=1MJI&?X9@-d`QITAmf6rJ-Jj3E z+C$2~z@g{8KXYprorweB?4#ieyg-yNvEi)Z$hM9V_bLy5@j^8hW3fM1W6o!1&`kFJ z{Ws|3cAblN8PRm;pWw150cuUq|Mvfiw3*3;$~SoF@b_EgJmY4fXvHi3{{9A#h^}yw zy&lfTlT&&?fwR24Y@>ZTR^`Ftv19ev)N~59f!;cwY$l%iC$qUin|q%mdLLXx8Ihag zFuv`V26&KKfJI)ved`Co2{nkJ%j)$`QG*wd8sq*P-SmtMDGv`3u+71b`_ArG9bMvsG1_e%p#ED?_$$_QkdLtET+9$Y9oe_BdW(1JC zr&TbLqWLT-ff+IH#!{ugso87fvTNfgnTRn=8;rGlFtB*}_%KSzgy{rq zun~GOcZ!hYbg}o47A8Sz5A#T3ua218>?i$?SOMZp^=w*|NO?v#Y+!nL_^pThVaZo3 zpa)n)9>4`~m{H%)aqj2fyZy4c1nH+w`BB#kT;ab!0rV|hU3QQjwhOaya5Q@Dj^d?s z21Dd?{^0Bji-I9|@QClhf&072NETS)kypS+!AF=q&O#cgnQw#xNB{{11-7he5?jbn z+tD=BX|Qqg%a3_^D#wc6-r~i!<2hxb@(UyPSz$;mRKkh}3j*$kSineX4o6pe4!#oE znYFEE#oCmg0rwf}vvSN0oRSfT(1U|OIvYt7lbO7}TW`o^q|-Aq$sc@HBn6h^J+Mu{ zo#-JX*4EaTzL4+~#*7U(%?2N0sotNZBL&O}avR{z z0IAZl8pykaJ|S2(r3)0bFifvD1A_I0{JCktqrLt87kz!|&Uy+h2_h~^*nl@}G73B0 zhHSy0_8?3yWo>)=JoM;P%j=6qlivbR%nAJrfOpQ~5)u-8PVt0!fSZH@A{hGwy{Vor z4?b*aw|)Okibu@Ymo}T*WRi z5ixhY@+Xn5W}>D}KRJa(qiNOW)ME+!Pk?FA#UvEvmyRKsL1{`Dc=gG$q2N>J+m4H~ z>s3-Cuim}$ziIZ<0Cr1)PV*#?rkQQJ$+t55`zEtV|6E*MB?hNq2lo6A(mXsoEXSk| zfdPCCAF8dDgvBNWHrl~*InhhZ-xPgg2d<>cf7XC>PFF1>VL%5PTD z(|cKSyllI+F_pTc)lj7*-TM7|-3xNo9@hE47N|LE$TIBgxTU0}4R6E%dU8}ZdkpLy zz=t>jW3CsGi6%N*8DLMK!!eNF-d?z)rrnxTucKLknrC$elCUCn&v0D7Uaz6=g@d3p z2jp!$Pr5FpMH;ss=tU*7wZIOq!+Ymu?4dq_f}-2r$}%$85Ha#Rc5M(8u|jE1$YPM% zXc^+$jx1SWpXpNz0__J6T%Pgc!|c`@x2ec zPc)VJgtbm;PLEds$?855(8RZN``36CsZ>8DE;ZrRe+RE3`ab}#QekwTXKsLvJd7w* zB7%(`G65Tmfl9ZVAW8M&DFN@wlkEIDj>6CuHjMOLLqj2GsFi<7PxDov^B)ii%pUr6 z;eKEyVE1MS+$3;_udM)@=X3&WQAS>#c<=e9Ic$j_iNOQx_?fCc)oa$j0f~QXME@;I z(^Y!xVF0Au?c2BQTzCMU`p5paw1S$5rUWDM)hl0!hJR~6pjWMMKCk~DNp~C+e6zJn zh=3pg=mid`+wVkQp5xYi-u`qF2R^*2Fb^Ad z@?c}b30rv#1af{(6t;qC1k{u0-h|&jc%l;Qfp&rMu&WO~4=7&eK>i`b%%O(Qy?f!> zt&r^{)KB|_>w%c}A(`#NCNMJ}VbMnn0Hi6>0TD!xu^TH^*f~R(H~lYc2IC((Ui?4l zcx`N876EJpgb=6@;s*g&SXfHdqmee4im>cqAR;a+4Uo_F{_%mW@&*4|vqZ1Ji3tPe zbdcV1bZtPE>F8*N_>Qu1FtE@BrH))#^qE?wItsnRF039@ud|3|u@w!flwM+Gy?EgQ z0^4Fc`PTr})c>=#*I3h{DD(SLOin!GO*E1HeeVAo=PyvWW}b_g8RPAhd(TTs1Ym;- zXy9`|`9>TbS)SwIe0YrURpu<`r*3gy1sNu)_(H;ek&}}sxB4I> zR^cvsgQ^&jumE;QoO;%7{raVX^87YzFmK(iMbM(UmjlIrFJt-NE7rR8jALy&aWpY5 ziL1W|!M#XUkNGN*Ap(#A(BFCO^hv{gmYxY20tKbRO-;&}j@s}Q5o$xnir}jaiwg(c z%k_E&_PlgxzS+6Eck2f=cs)gberFiKa+H$wBgL2x%um51eT%8uX;64b(-y#M)Kx;& zg-3Z=q{?Q8EzY4^Sc&R;p^|Pvsn@SmJ`i^#ZdN5ndH(ZLB{+JdZc6Ls4H4{uvUy^jo zKTXY}9rP{riYR!koS)*p3+qM&^zck=YedlNbWQLfQYj9LGqCD(?4 zdAq>@zfWt|vpGBh$=$Vd&inVDUGz$2!~8-jOsy3j{hL6?h6nw^g(-m#RDd|a7m|Bn zBRCxKVi$#f>OLyYuifmb67-+`h&M%%hGkA9fur^6v?mQk7 z?ez_`VM^1|Ua9<^+t{rY=4b79ahQ+MJ;}IEtv>G9!yontIqaVm!h7@W(WWs*2iX-i z>i1eXlo!op^yOwoihbaaTO1qBd81&yyijyVlWlG1R59?tG3>Hyf(}^(j8lDbW{H9P zk9?RsZ&Y7QdRCcL8W8`G794u}g0OQk*7WUz?Ci-x?;r~?tr{wtENPh}rgO`qL|nsN zMsdIA!v7gGEBh2+F$$V&|1$w)K6CN&X4ZKJ_Md*r6}Tuk&@B!8qJ4(=65w8Ze&m(6 z8mizr_ldD>R!-PS$@QDxrIl|8-9>mjuA3>2enhN-!slK^Q3k^I2%eX#7pe`Soh|Ur zEgtmemcUWweBp9XVb2t2iJ}jbmr7UWg0i@1QDne^fAIK+!S?j+c-_!?fnz#>l7<83 z7X1^27nwmvxFY<4v|l+ht)y}vJxi{=7}zfseOohIG z&{9md;h!F8lq-rmR<3O=CrWEIJQJ2-N-l-vX}GfT&?ndK6BY4dsGMIuo;|a0TGZsZ zHg6H1GV&{iauH=6CRKNodu!~whq@MirEIO3awRQ|#uj;2pW0gN5Zf*8(++w3GNx=% zJlKh2WFJscdz#cZK!EY=)tT!45iFW|=6@H$r0(4oa?cZ9Xl*R-_}trj>ffqBApNSQ zBb^XSa1#lhxzt_al#-uG+%0v~+p#5qVH`w5)>r?r`5Ns<$$=3rYNfXBa%K{$WV-FN z9(KdZ%-?UcO7#B>uLHiPNSH9;MBTQg(sG-Q`=^6iM2DW_ARSS$XLU|0!9k=_^bfDK z?HyW@iOBtA=_^UiJuga0;-$lTqqMZ_t3Q%so;CGXrdE&83!H9#{|TQPd;~OMpQ!T# zeG6R)Ha0dusfbh6PD{bdOetK^4R==P6P8)tc4cQBoGV(!Pn-|^WWc0YJ;ZX9D`Gkv z=L;Q2m?;X)4`InUVOOcB-+hBC@`&M!v<9)?qkijuMg~>6-t7;r6X|M-E}yx_kz-LS z*s3=7B1IyUZ@?~ZW>dra!GpRP%8Nke+_+G?OlVBZ_t@Gaw00lw#$!QmojW0j2pDCk zEAj3fHlkVO&KD9ID$iF-FQ&<0y>32QF?6&NLt6AW-r8EIdpckId-W6kx_$*!WeW|l zRx2@ji~TC8_>!pEMrliJd55p%sD<*(}i1s^5$3*gd`1eu5$X#&GHN2n&)^6%_Kx zHLUs{K>i|d@8Dg9EeLZ&sN(&p?Cj|cm5eW$t;^U>w>p(K(ynmt#p*ubnL3HE5E*{s z{CKo4VtO_;fN(G1sYk?DZ68x!rp}^Ww!FULm%r_gP7)ZF(}-(hrzgoLO(e29&3Q=hBXib@rD4Quo$;>lHH+|}!05J_MO z(4(-vZ2f!nqLQfRz`(}k@J;}`S;Xf(X~Dy&Ng-e$rlFNHnGWR6RjY}wtbzX1snV>u z#rXAjg~8WPwgbvn5T|J;6TGKbf@)7BHGR2Y=Aqw!l7o;)k>I39yt6%=r{|;&!;emv z24Kw4yNX}z3t%!BiTUp;QU;5i&5ut{s(7*)KbQ<2izT?k_X?i_s;eu@tgZNSQ_%Et zjNs_>`o?m(qSVZyrIT*5`mhDF#o?TuT<1sGZRTRudDi(ktcwYXUA>dv`y$nOoesY& z7W9;uFrbOXgZH?f0{181VF+4t8{YNx1-ccx)4G7GXyKCK0tNNZy$^?ynJd1(ndb`V z+lq=Xuv)1B4-T}E@YL)Le{^SPJ6p}HUPviWSeboOXXOdCiVXtPW2327a~QHsT`(4T zGWT7@q|j3@6yTS3v(*9G3`UNu+%x$a_4exb;|>OoUNJ7=S}j@5Iw)a#)cq|O-#bgr1d*BV6r}A<<{b`dkx12_8#uD#65)NT!%~^exqys z86qB=Y68yRevSEY1YPSb&*Pj{X)fZ+7IxnGi0Lu+qE84P*|xqCRBV#I*`TJQO6%K}Q1r4cH~V6-w@&}1O2Y-7i>3~>uc)qPhh~Z82%iapXuwe~)mNo1 z|0I4MI#wX54v~!C0n(f)!A90KDR=0=d5cr_G}P{NgNRa28Ye;>-;(AUP?_$%9tGVg zpD$b+U^?jI1nL_P)a>0q)p?UY_xe8~#rTjt{+gOSEO5b4cxb3NA=k&Q+AZVO zXj)Ob#>@P<)w|7~`>Rzvd?F8>v*!MZ1t{Nx&@Ms|-51zx>Dc!#_u+L-qqklYFEi6D0Z_*~a zociZlkDa@9j;XD7V1~B*qP?Jr-8*%zvag-#FJ;1ZxqcGMhNO|0yr!n~J>Vz*@iYN( za1w93vvZ{TNBRWGNNFob0|*hIJHxh*ciYp(ma?0>q!;-J*!0BrK4*7aR^A)D?q^c< z@i_1#?mYyuPiqKFUkkrFcP)svit)cA-u7}V822O`U{yyZysr|5RO`vg185kL28Y=t zx%;@#rmR+#VXNmgiB$JL1_Az&y2U(8b7u_P!#mP@~-YXXA z+nVn3_pacCGmM%q%omqETUqkDC|}|XVOr;gsBtdUW4H3P%>~K^e_SIyCw$M!D$$J! z{+Hsndrl|G+QglE>k(-VGv2@t@f?0)FZ)tNTDga$3vL+rmi@9@xw9n@H_5G`Yqb=g zwbN%za|QTMf7pSjP#8L5jh^2PFKp#juYcILSZ%=l&FmGIt_P+4z(ycG)@X_5}+tf)49nXMYu^t&Hc_G7fNu&8E~DulMvB;L<5M z^r%#&T<~G?+_`t1!k{JgXP?-w;<0P<@APDb_V-isPgX+VilNOdMSWC4cknLe9NqU% z^X5wTS(RJVuxa|bQW$1R+&i{2=)v~m!S14_5wR6DJ6>e!%zW!B`(7-XpLzaGp@};C zZaJsr#>H(X7CEV@u3HO<#@|a19#5%Am_7dUOkHn9#xl1*y=POCrro%ayib6vuPQd3 zgrVFWrBR`eyX#|LPG|kF)016p*{7?U8?hnF2#N(NOHFFy#_oLSgZ^!#RH1 zipj;r(5uh9gJVQ*9p683{e9-}?$dMJ*H=8Dk#U2O`Q}9xQ`d2)0}2Rn7n@odDQ~`f zsMFMb^J1Khl_p->Ryu!SadG{xvE_~NBBxpTbdF&IA>5ksp0oJRgu0*S{AsTL&I4bU zuvQ?DGxeibhv4Jnoxhtzxu@@T!Bnt5qpxn`ga*sqy&SfUzoI4XWELGq8rx+BNsjJg z_3}Evr*Di%*!8#$MXomu;`%uofBG9wrd#KV+3lmAoLTmpn79xX!7q4lu}4d8DB6Ml zmFT3Z7p*vM^!xQ5a^QI!vK`fo8-z7ID=U*E)0Jz>9182=!&%L+p5%?Oc$hjE_>#H_ zYVgTr`CWE-DtLhc77EdwJ42IpJK{SARkAxje%B%w%jI+ksvOK#jfwlT>6M9{!p%S+ zp%sW@sHSJ{i6@jI5-P_7X{9N=liJ4l628#yc}93cvUl9yfAS%C9c=FZc#p7kJs2vc z@Jby`>h4aEu#k4*>k=rty0h&3XXjo211$zlc%pVLxIk#-1faAkfi+q1&@Q$@VYvJE z@7nzLg02>Byw51m-bA#^;l!z!h0tFcI zALi@V32^-kAC?3FylJVp`iw!?spaFoH<+H&730U8zSqTaGBU1hS1l}9$jHh*@yU$zp?L(S&{jkMP}c&>M!tEjXX z&dwT5(}*|}z0xOsn4zaFR`J-f@cmM;c4!0@d)D%ks;oCrYD>?BsRc}4iD1=pZ+e_0 z2@!ISJQUpQ?l|0D*bb?Wxj8<*>N4y3b3eYm<==V-WdCm~f&VG{d6Ch8$$Jm1S24M8 z+JLF4DVrK_6rTaG3y3@%nNMZsVrQ3ue>({?04N~Xbs0!XI@)yv>P^(cad&?x74cNt z8?Yy|&Zqc$#OBkw_#h9kh3 zh~kxQmAl{#d_Sjv`uv3ZeDGDkYqWO2Ptm+j$dIkF$}~uk|Ee?t8|5G!EK<_%_12$B z2ywACHa89PH4S)a(F!#v0Ql)sQehzlzzKkCLCuw*Is=so9Tt}RgTWFFC>J~puOO)5 zfR<7La7nk~>8;%UY1FiYvJ#?VV(Z{F*4Ec!A%AhP@PmK5DOxB7Ojr4GD|BAdPA<%L zn*iBb);}IYN~(JLSompqI|U0VX{XA6J-Gja5!GdKs2dCwCg50dX;S}(zepIJq2~`B zgoO);zqH2j*SVl9seje zwE;B?Mum$)?c@Q$wTQC&uA^j)!I`_0$FcIl`Zd{w902QKtX^H|TVzW`{1l5PrU*ROZ ze_$XBC#TY*zJ5)7=CuZLWk~>v!2tni{C=iYlupqMY6Ah$^%`tDs6>H9NJ!Ya>fqpj z4pIh&aOlr3svD`7Z+8w3UWN1Fjb9$5`8OJoLd;73^hsv(B}n zP!)i*Lv%9?>M8*m2Gr2Sz8u|5(-`RusNNw5w2vaYPN@n27-)qUm^o;D9H3FAz6=eSwnp!|6PN_G|RF(AB9$}8l&F1=sf1@ZmgfBtM4%}vvE zoxuqdMgeFcA8hY@yGbd-l3`gDl~)^+OG0`p8Ni>Qjjk9Svn~spzOSmvT&R^T=yOt; zr{b4nhFBgo3)vGON>rP^YE=Y)JV} zO_SfQ%CJIc0Bc9XyU(Hstre1=OrQbF8xRy1VYEO@k_N9Gxd25mz=^!3At0zMzJCeW z0FyA<-lA9l;~1Ykd$zWCihx^iJNnz1CaThw^zPs zcK&*lt&U0r+=iHdH~vLApzv%#!z}>Wh>LATL%=D6m46bLka{Ru^a4y# zr>9RbUzvFpYH`5-|60Pqc&`H@^Amj(N65xc4E{ z1A*fM!2G{=(D|$w^BhodI9?!)|FKZE$=)93JOKTbfVz=L832L_U8wSzc2o_>5A+NN z%Fsh?n^p(EI^UxsPfd18c=jj_2J@!zE*Hl6)ta(YaQi~RZX~6lxeB0C?o%P@4LG!L z7GV(=*MoWx3g`=}b?5@r5(4u8Di6s(UX70Np;V&pfO6%C2DmQtTt&-ZeO6B7MA@RD z;bFBL8#^KiqE3a;il%{-Q8{SroPM7(h^IJpdGw3Q|?0rZ5^K zs2y@m-TTQ}tRPzGCGKBA`Ts$-`#=3_eA;ka4i1itzC0*A2YfKpjB#oh+k^pB2UU~$ zm**}kMCxok%~rpT-uVxp_n~qIM0Hm{T?~q}xZH=u#DHgX1u#0uI+m+)4f6y5`|S$v z4sRzDlq`m$BOca1Z0^+^0I{tD!g$V2#A$!+VFlu#A0HpTwncxU`y22>l>H8YK@Jh% zFs{|Kpsv1t`c2I$T37`j0l4t=8ygqkjiDk3iv0xGPdC+B6rpv?@-*3D%$I*!T9Ae|m2HbFB8*bu`UL|G)rWYhgK)hbbsoK%X(yyr=Ugf1F z^gdqHp9bo$76s7I`g^B4wx@w`c<2H@KP6bn*N~*tR2*2Z^9HRDUpUZ@a^Z+C79P>W zkYz3o0+|guW6(Xp?c22wcDk33fiy%(I3TKE`NGSs+&GI;3nNb5_eufrKtV3R6(?c5 zpr^m!+!6q-e&4=*v*D|OQzYo&9iTOyu8xK4A5BrAb@gBntpVe_yWodbHG;!6%A)Jo zWaQ*fl?NL!x3KsIb$&kfFQWL|B!Hx?mS{lz^!c;V&p4sd<>KM;5J=!Mgk3$~*ILs{41}FBDQjGGs{V zX;7J-kjx=Xl2Ya|$t*+WIf@c`ehnJTlBckcC}e7&C_{#2CX$dbGKKU0^n3Q%`?vQw z=h}N;zjK{G&L3S@t=8~e-}M>p&wan&uNx3Ct74GYw%}}ot^PR*=)$2Yr%!X?aseTv zW&pPWXPD5`;{t?-8#Zhpic=!p!>rp`>Q-X76LF@Wp8_tH1FCGInG##5Z~$?Kl~60!U<6;97vBkzdAv9g2Y$?5fOxfq@eci2-v*UY?nmS#0+^U*zT; zI}Xi(Z}#ELx0W2DajbfRa5TrZVH;izWdr$tlumX?mqLMy*kcK}$A}!u%B$FFL`IBC z{0gc|e#yLv1;oR*8~VJ|E+XfKygi0IADaBDVH=FEhvB?gu!>K`>rPW{9U>k>;zJW$->DUaR6S4tJq=H5ny19yuqA; zCnOn6FcAD!MZEc#{tyr8;Qngi86o_at3fH^40p3 zm(&n7cBfaT^-L~7M2vq5Z`q<|U`-qf_$ComRYzc0So^q! zxEzQ`I*6M76!0%|Jb_Y3{ehvG6Sy2;LE?d*Tf$NX*t-0s#8UaOY^znSXM#85)b&1n;=q+OyRBN zQqT}3eWq-XQy|~)$Okfmi4B8Nw*0ep*}UoTX-FMM^A<7k1Yq4{q(HH}(~_oxO*YgD5e1~-WytmX%UgK>}|83RBP>ujWew$Y_N zru(I-G&kRb#SWpvO{L8Rl4Vp020?SWvfA3?bIWAd!G^LKyWTZucE~1s_jc9aadael zJbU)+rY&0-5W`k}dKKpKL7GPSdXL(XZNXq@WJK1vZ&<{&-Q{R?dLQ^AHRhHvvl%D~ z3JN`aeJQ(V#>Z8xtp^gj8ocEd&31*BeXmuH#_jKc;8r%^r!1UZz$?NnxZl=R1QIn9 zr_$1;si~>ZPA9!)_xPov>-Qk0H^f{ZE)c{cKy7|V%%{ZwG`aG8>7op%f-c1K0zBoA zo5{cge>l79hU(QL3pVODNtpI+Y%ChY zN%#fnC`6Y`#wD^lYF)auBQi2FYwb)B!C+6&S+yD#vQ~NmuwH%Rv?!A)YsQ=J`k3kT%r?2(qN(+ZaVCm`IANXFg9X3!bU=< zvba6I+K%s$+W>kF0b0Eh&3hz|MtEG75VR443ADk3uEoS`1;`tj-0;R-gw~*lF+`1A9J& z2~}VN@Hv*l@h&QmN#nv@A8eWvgJ}u9@K;%R)mvTBGbgaNT)T5 z5z8`S?~o;VwLFu2I&Kb*U?e5;syA)^+YwjSCv98T@TL7X@PUt4v>@xj4jviQ_ErY} zaPaR>aN<$z0hKEVsVPJ%I2}FEctt?oooWvn8;jyl zH1A2h>z8?wn|lIIJlt7Wl-2~+tE)fYd$vjf6ZyppgRh$tREm?3PNZdQWTd^qp99N^1kzhtp z-UC^Wav+^x7|~5R`&ZIw+pnb4-Q~@DIO^HoUHEq4z6cVk$w^4fV@UnCx#I6emiHe8 zg72eE{+a6XJrd3LBt=_zs!u-LTt+)O@`QTm1^t%Id>2KTYI~l(yT-086f1K#r@v=I z-qFG{6D-%lZJ4j+otsT&_@cW@wQ{ap-(I)q!7lDBWO*-}8)l2!eb8+gp_l%^t-1Qg zY4ShFsM;{!PSt#KQ$;9tuMx+FkM3;O8q4url7*Sar(YJo6?DeT0l;qF@!Vy6d^{*EZ9613CmG5PAxn#auTG?15Vy6 ze_;OiVK>?#CIlef=5YjS`WzzxxXS>_rNjQ~dK4hchz zUx0!^&gK@{>V8&k%${D%p0AsCAH923Bv?uIeE#O&AMBE7{OotD=g!9m<+f5r6Gh@$ zy{)B7t4kMpf8+bm>{Ob-{Vl>Um&>C_CA0PFxQfuNom-jkfzy_WyY5lU#5?XQa8M*{ z|2FTb9-FlevvGP~>fmD&3iWIHv*1x{F?>|_UYVp9xBMsh*S?T;q zY!ib`VxEq}hn}CB+2idh>*~e!jipzivb0ZMTESt&&&P>UGoxsCPxKCV#;6)!jdACg ztMS&KXO>PW_!U#zgB!g4g{5STTR&!+jCM5t;rwXSdQ$20+dsarWWWfayf}CqwTo<@ zX*-3P4mZ1)-vkHJshr9JWz_gy<=*491;ZPLw6km)Wg64f2BsSxMpnj;)%3n-ZOb%K zY#*Xpx(jOz1w2f5+2EjTsvbIkoxf{{*UsqynD;> zI4d856t?{1jJ+*2QI^(p!q;U)nIlKy!q-f^c|{wB_Gx9}ScoAG*7w@vCuSz2OdW-4 z!nTt46i~KQ&CcR^KfW=$D}I>C*K<*l>ejh8)b~~3)3iy$O8rYw-%n^=V7zBna=rHb zCElR{8wtn3R!_z5SM|zksm=1V7Dl0C+bg0cUb^ToN$~|wlyosqO<(@8)9h5^v!rr| z^3Wj3;-Qk3rIPjjevhK!w(P%3L9j1G@D)h1(@~Ifr%;G{bl&enD&L9$!J&vT%{mqH z(Q9=9o?-*E=WPaYS)l^$F1yF@t$?`#qePXDXWe^8Oqd|0ex~I{$jkyxC3J)qUJ^ z?R5P1HKpQmKVES)L~j>XWN7l1ZyLDIE>XdmyG6~n$)v$yFz@)3!r(q_-LOI7fW`Q; zvSZzg{^uB33d^I~Rko`}(HvUl@;-8R=LAUP>5o_6w=`xvWkxspRVXYcSh!d#KvPsu z)`%8b@k5IQ%Z`TZS4*hnr&EGozv~-6mdVu=>8O84;t|by>NwLb=UZ&-dn|TW%ti^i z2F~zAg>2OcIxKlf-@t%a*ZQq9bt3}di{3KXRntHh=OHgAFVJ0BH*+)cyS?4LLJwV& zud)uyG=p!lx~!C9r$n|X?Ms&TA8y-I;ce?&Qokm$F0n9&RWOboADZ3rO}5>9zm7eB zg>SpBqrCJPg{6s31A#{l<;f2wo(p&QmU#&6(}~tH(t8`bulw@DW*I4?Z^d@u_XWZk z(tr?ZyL|OKU%xap-q}sltj&R*o<7*}MV;O7P`gaT>DK(Pn3#{;%)Px7Wt<)ma!u7W`_3ADXA-Wbm`%X*pPQQaRqx+E%KgP97&9J zP5ol*`Bhdqb{~zUr_y<@u<-D{!F{H!P1|cmlnok^9oi+E8T`)~n@q^rX~uCkjy+#R z_fu+_ecOs_C$vJNVr_Nin@=a~Gu5%q-|1*OCoJm|FnW*gzVrD>^;C_6iC(Rho2AOk z28Ea^hcjj5GXw=?jnBPoI>Bi%Y@-&Qz;n^B*c^Y;FSEHxlYZsI^QV(T<>?gZ&iU6l zaT3eaX=f^zuKxVD^8*2{O7XH|zIrcT8GtyM9-$to(O_=?OW3baV%v?ui(S={mn38) zR$zEu4ve$au1Qunm{)H9lnpD<{u76{N>NNj=eds^6P3ujprS0j0&ByTx9OaDf?dpl z!y!XDO?P-@`ONh8CS?}ehvA4Sb;NG>p*Zdb4g;tvk2D7VY)JLRtUBTj*dQnlgJ zQ**Oe89L#m5@ln@+2q0`K~9E)DM=r$P6Tw>*bS@W8=cylU+$S9R8`Ytam1X)z?anL z*CTl(WN!W)M*e)PoQv}dGka^jLuS&={K1w%SPVi}PEpG(Z)$jt)pw=bq=iO(-~a9| zHRq0;g@S^H?R6PPg|$DO>n^k35#$J6tdR4H{NYka|NX$xB39EB{~r0~ne3kn9aEfG=!SA!>YzB|k{r6Q(H&YViMm%gN zQcsOcnn#`w(Fesm+c@O0)uYQW{aS&ge2pVrq3Utvri@0DhK}iZ)OJJ~BInJ|3+$!a zx^1---zL$`%)VL)dP6y@=~(L=OXD^UWf@x(&c7!k%TEevx!kgb>ebt8WJILbN;&b> zM@<%M;1{m)VCYWDUvj148SPQLOo6Az6V z<~us@BYLJjk=hy65E;iTG;O&o%R95U8b5FSef+b&QG zfJdSfiiL<*+zzoG_YYZl<;G9vj>h1?+i?dqG-*?>t-G5kly%z|EM{oi^fLvY`&=z` zWgPG`E6>Rd@TG1rH8tfO8;HDfWl2_BH%;I*=k|q9UdNc0Y!%p)g~S`5dhB45uH2B+ z_pQS&_f@-2NP5M5S!K#_SvrTZR??*?e%ABHS$FccX}6!b)U?f{#xBxr&ST&aFmpl=!Ha^V#d;-((!` zad0MUF3RRBy0=9)Z2h0B(z9%aCRRu0)x0S$ek4@LF;n)Jn2vv)QV(gR7ZPXNt$C+q zJ+1dIGMW~vh$U;T6#KVQHBdcoClTTAK8b;{*b<@CGL66=f zJ2Vt5M%i{*az$1J@P8dnnM$)Xj#J*VP`cOhOd^{qLuE2Gk6}7f@JOWlm#BcsOia$U zP0A})%x-afD!%q+Q+c&o%GAh`5H&V#;r+zDxU)%3m1nA#R9&dnq zgi!ctpTFc@xt9K=FJH;JC2(9*^tRJNyu4ptt1|?(|5{c=&yX!3HvmOYo?Vib&NhDI zdkWo!IKfM^i_RWkUubxHZQfKS6yOoOXTXp|Ol)S#E6bvFB-8ATZm}8ku`xMYW*W3= zwY6vo>_}I0*YV8?F{ItiE2s_QmEwzZrHtJhqj_G6OJzmiZR0Gb!uDeV?lAFRp(b6; zoLWmR<3KmVR}bwAny2MdJ5Pl4DKe&%<`+J$EGRWuIVN&+&yr(>`s}>d%UZFXK4}(l zQ7!HK23m0Oz~RvJ{56f+qpYH?PfhV&QcTlir<^O2wQxm*;4WX2Ex7R5^6vLy!Gqza z_U6s-Etgc3?W6c!e(JU$xP8OQgXuRHn;R@oWi4G;dLJs?cA~TH^0hmwv}UNW7oDze z*9=kQa1Dq|5oK6;6E6tveVh>6-R?5k_d>bQC(salFC_>WgC4Zo_)DrXK8F1KA~i5+?`A>y#^#f53C9_xXUq?pv1vT_ zIr+Bfl>nvCvE9Xe;OEU&bA<(3pa|2N`^CG0;@Ech=e2)286f&m<@f|^SaAIN!s}_p zmvbu+>i@ybvBFHKHlVWlbSFj3+emDhGk`&{<#WpqXQguaI_uhku0MQY7gd${Ul}B5 zvbh}MG*@%aq5I-h__k^7$Lv`RA5F1^k|}i67q%4Uf$>ZbN9~!shTmc4LQg$p-3|wM{htJLSAp zi(3a4{8c!Yyl-8>=S8IRF)xFVvYXR6OTYR3^U1-MD<(?{#aZD2`Qc+0E`4&W;IiM6MsG?YuSq zy2rVP>&}gcyH6NC$E|ao))(hs71!c7Q((k{ok_A=Rn(1PW-W zH}U{0uxrVdj)Rp)OZa1UFJ4yf{MQ@(vpdgqY%3!t05X+<`9pmZ`PHz=B2)xXz(LiH zuyuCoBBe0;GIMEvb;9B|Yp?miv?5V?sLmkPLN}+2DAkZeI@)Ckq8@PSF|?N*!AQYP zr$VI&2mUw!r%!Dgpy(mm%)Dx36>cGsQH;VF$r2+i3GD?N1ws;aFR$lMqIA;_SXreS zl0?aU-g?S?^tmpe%tSbE z5Q0L5mYokZWutWROc{JYjCoyZY>J@-Tc<+f()0E%m@N92WZ*E8fq$S%O2KPO-r7OK z8CURO{TX-nyvb!!O~r>(C`jl7^M+cq1Yaacfr86oXhp8j`Fe|gsh3U19#*U15|OlO5>t~yQxi(vfRreQN!LF^ zc=o^1qy5t#j`RxhM=LD&zImnX*5kVrEs6om1mP_4yyV%(Jb8f=KmJ}fcmegLdiVK*dVLjF@6Uuxpo9p_M^CEyHl`Ad- zKF{!!uI`V5Sbzx}Ub(K6L4_K8)qL6{+SPFY*yB3KC#;@jg@m6+QD?ua9)fffN>)?A zsZbLE9v}@SwxK zJc4OF|A2$)q-4T~f{X8B+=l5fPqWz9wl@&R`Dr&e>^~bWzOlqK&BH)l>f&Ex_tANS zs#EHhXJ0mc)=IxZ#ZH9u_ zHhN>z5euUZ&i2i0TEdvKeAC<_Q5Nhd`H3IQwZ*~1A91&ZEieNB)9z{r%~)XV4o0$)$>71|GX<q?>PJGUgxdTH zCMOa0gkw{H_J9DdQ0KS1z3;bDU^-tvv8sLz@h$jOgoyYT3k!v|Eu`=*UCsY=VXV~b z#dNa~lmWZBBmT_Fig7iDA|1-aQ0O~U!>JQs%H*Evn}o1999kSAQ?fEJv&=zf9T#-O zTm_*OLGqo{(K-F7K%Cn4#GWQ*k!OfcL}wvg+oOjM6IBj7C}6#WfQfGfA;aU7GpW$% zJwad2954?12Y@*X$Dcp->)%N7uKCp{5Y4mLvHao^uHTomyj*pQRSO#fN}D{$escpT zMuZK(JBkR^#P{`6`S|${;G!Y#?7UqqD~2k(YI@Vv?*p}~)0=jZrhUM>r1JBm)VGWz zcGL7x?u=klnl0=18BMk}iY2 zL5_9nF5cRucLr5=eHf5Yz)uy3!4}_#S|rK(lbtgth&1y%{C5RmYxJtX*;hi&-#_eS zxLztsl_-m4=3zT)C&kCLvp^O?2O<_S^`AaHyReXMQnRuZT0<0%ixla156|Le_tqt; zd42xlAnN57*97Ht3m4MYqDTyTXHdP3I`W;kgdDKqOo!h|h=m9bD=|?Bu^D>09KShN zitDpfaY$ePQ78Vlt6Z#cO?Gd3-$alnb$AcNQG(EBlmLk+sFXwk^&x#!;3g~`1Tlz( zfk@p6^-*E&^9@p0gwB)VU2!nI2|p82F;j?YSpZnDDi@UKUm_%7U$>54TwI*IH)(42 zm$-EL&UNr1q2u`kS&*LUk>>oVFuE_eTPxA`6Wc$!#Hpa9zJ@WR&5HFNAV5m9o1ixu zN6pFDz)6A^uVzZLuQY|I&~WQoL$uH%1`MoCP&tA#+lP0p1VM5&PUlBLP?uYb6=A`m~00e9SJ?0M=nl_ka01_(@lo1`+=cZr4x?$YaKXwMh#UsBFak=_D`@IvTiJh!$CLyCs#&xMRL}MPzK|g%B%C|QX?4(9 z%)p`T>7&{lOe?bgbK%m6Vq8GMOS^xl$gSx%2K4Go-E;v4EGj7(r2&gRL~=<0olqbp zo#~HZEWHv!dgzhzDTu|t=*9nUjUWHnUG4vV1?~UqD_LfC?X&tUDDEx(o)i0($Yx0 zwf(;P-Fx4?W4tl$c;o%`92T2%_WA8!tTor1a}lPhEJsX0OMpV5h!y0e)ln#nLKF(~ z?gd=9vUdAcKKz62BB`Kp0e<;hFb##T@g3zKyP!~n6v!WppA&v3a8bxrM%(qFgSo4R zv9lS<&e+w_*1^@*>It*EnX`+PgZ(XTK5jlvW=mIBM-d*L|Neb$2WJZ&GIQ$7C=@eF zLHeGC=ZCckZ@o8WXA)<_osO3lm|x1XeliYCc}Bl$B5V1Pg6W&5e3!s|ef=?#zgY!? zWx6)LWw&(?v^4FDLX>Wl4k*93U*rp;RO{SHHGP30b2oUd{=>!;Q!Dn(dgi8a7vh4Y zxtYLs-ix)H2La)*9l8O)(ZqtWx+fuZtOX(Mpc)6<>fjj79t3qwOM%+Jra>*JC~ zt5RINxWL_Dzb47{Z=LMbKlihWZi>8g^{Twjg8rmgVMjficlGDW_TtBI-gjJEs@JyX zYXZrrmT!1GmXv>pc|y`SiqStho|$h{@q`7d=?9NZ$SBst_R7GDu>Gcc_iXaHHp`qdGwQJE;Z%D|^$=aGoC0^-Q zxw{?vHDdZki=|}rlg>t3pYoGP#*SdxJsk^Fqm^{niT0_kTzZm&r7hGRv`mc-*20Rm zabQ*=c!^RdgwlGYE9~}W^XB#BUg~FO$Xhu$&Gm(ZxtzV zA{l3otgI~Ua(n&ohfAe<-`dBo*HiHrylM%j#v-#3xuj*JvRq2(nunuzf5zNKIDjKn z;7S&g0K+v)3g^6ktU)b~uofUK}Mz zWEuYH9ZCHLOLcnR?^%VaURe45S!jy=M|V8kCKHu7Oxa@7(;vNQk-yhA{uG(C%^vTq zq5@o&x{aHHiAknHZjf`-FLoz|pZ+~aD=E2TT}sPSvr9q6t#{Sqv1{yG%=~mm8>^>t zs&rpS>fhan>uM}WKREpBDT$fFRUlKN@Vq>val^>Q)3!|8al@MO_+Tx#W7vjIs02^! zjpDdq3C}U^;cfh~p<08i3a4lry}G->hSW>1-udlsOc0_D|4yAp`u^QFKHgpNTuSnw z|3<0>vu`(Ep|(Mh8K^^8-r8Xoq{HN&Y|b(~>K$OHqm~$X=@)}7ag?{mLGTFr-g8Xz zaO@z-XBUGgD>NTP**tO^GO}4;Rbd}|&2~3k9sSkAn2;P3B& zDJXJTQnsD=MtME%$IO@qdm}cD8&{d+`l8`}*<}A#g?PizRz<7P3&pGeObQJ}oFCb% z7^nt6g9`Q3a$H>-!_@~=hhARw53}m=?QoQ-A4#TTbH-rDNn+sN^(I?rYYPU=ReU_CyteQ9u-f4TZ|)H{DN`v`f>!v)i;WsS1-M0OPRF1C&)v2 zef8DVyzI7>MuK1OZ{5pMt)CjcXne1B=ym@ft!k9Ewf*^Ws<;DzleP&A_A{KL2iv(QPaO9KQ{8S6k*-pmB7- z-6=Z#G)ga%*fvERqv-m27c4U~Gw1keSIApeDuhpFn*w>=9YafHhI`d*gz((IWl|(B zSgUF0%gH2>6J$Mou0VAcGnd6826dI~DvG9)qMPc`tcEfBrHbv^U`@RyldL`y751$H zsR=SVDwKo+X=b1gUy6r^yGqP`lK}GXpPe#Wrgn9t>KS-?Vc{DVkH6_)?I+cl6AQDm zC8X1F-9o8mUwZ9l|4rd%k8#G~g3}z?E$x8&>7^AKI^HJwD8I?;*vDQ^94$!?y(WY- zcDD~2=8b*HnHcBmH?iHCk!V+6U$5(@ZwgODZKe{xusidcgC!$}KSed8&B_iZUt=gW zx7Cp~b<4$w-*5bfI~OYZcNc7xw7bK*!4XVmhO$pQ$uQ}NIs5233+3y#&YsCL<(ubp zg)#aa&XCk_-^8VfuGb-`L?kQQ-Hxl--H~FAt@zi4CYKjZH2_O4-K*mtti&tv$jF{(aKRZ zOZhXdYB<)vl~4hHXt5$G-=uwE$kKGc?WWWceEjw){;!xc8Ba&S!>2jZ3#rX1uNftS z-?J$+PSU7<<>QsV4HIR}!P!eq8o4%Q`_$e~xc1IdoEFLlJ>hY~G>+Ms6boKy`eTBM z^(qpiG#CLNIw%k4L#ci6zPv2Kgj6mIrTynVqSdxdUJ-y z(>#DYdII*84$&K!uCYj~=)H&6x()MM_h_+3`?0KW9V@F@X!eBEw8J9gN|pz^E*gAP zMBOb+7G$I!VJePo65SQ1zn7;el`}}-y?Kv`i79i?y0T-LENGr&U7il|dUR|otH)(& zeRU~nQ+I>QFMsR>NIs{iugu}k%oTL={EIPp{RQ6NN~MI#i_wuu@Tku3Tk}e?MOdwW zKd8{lmXkww^*nrb=d3m%X@K;;IxB;wjTZaQqF4-B>90S!c-pQA<5YN-5)AIiY7AUf z)9b6f;WicumLoG$aIA*dfpU5- zRT2VECmnhB(7PmQ5%-(*q4G3|#XChT_X9i3_%p4SWUMEz8zqR4g+H`!J51>9zSc_| z>biD`quzaeH0VmQ4~2R0Zf+-4D%SV>HeVb*(8H(_?el71?$b56j*p9*r?8Tx9}i3 z%jCfP=(oiT$a5-0NW&g1KQL&nb!o}Wzcn(^<$3gtamA)Cq#%T%m}0`H)+~%+;pVv zFT^xaqmjVqPuJ8Lgimf%<|SHX`WPqb%w?bU%*CB&b^t^BuMuWo zr_5D`ntL?{Q9DaLqOOSw3JMjT+vq=k{uuGn6k7FjC?$yurv7FAl+S)g>;m4$8fyQ9 zOF^$UKDy_*=;=OFyK(TOcSn@>Jo#KET6yaGBXt{kS|dZ>J<8eHSzfE3Y^eI}xi$uh z+l-8io`2WN`J85@rRC)2&i-!7DJdz1l;4*wVdJM5@qdr%`f*$AHk!WH-Sk}_CjDU% zQJUsP%Lk9T%@g!@ju4VG6?ysk$r|tI_wQehkL#i83yqr`XPfcm)*J27#5nBeVL-%l zm0S`fo=aA`UPn)jN1CtBWE4`Y{95UyN(c$T!^Fg-uT(&zGYAQHTaIN|^YU(bzZs(0 zS}ThU^|Z&+HB+ALMTJL1v_JwQDjssFheCdxV&7D(Q_ID8ykB%LU7V?o;x(ubldq@^~3pn z+S@n!d0F$1u3G92F?F{w_>Vv98aYbR(;uGibupUqJz$WHpv8$(G}Cuq`^87ZbZ1f_C-LS2OzDZS`lXW^U5O6{IebqBe$?1T2)~xX;=F#HpyMr!d=QFZ7q~}qzzME)uB!6qM7}hcyE)FqqZEE41Znu877Z*cEIO3#WUj?cXPR` zVz}mQo*R3^`kj26a){DmmG$6tv7OFap0iak_)ktxMNM|>D6%dqL4C>mJ2CccKc{4*<`fD3W~6K zHzR3jflg(T5OXuylbZKT!FwD=?aA)#s&U$I;O@WWyeBPKK6=jLU!tO#T`i~!mQ20t=_v&GaX8h^?|3C6G%W0y zde)`IW;D5)ceNPT(Bhw)HTWwIOwP{yvyPSl6H2M7^-XCeta~OVCPj8*0-ryB7RtPi zMi1W<6|Qd4Et$8z;b~YzMN9u2vg7MzG) zs*qpjesg?RwC0ZKZ(fuOv>Pkq*QFvO3rLj;L9?-aZrN*%pbN>-DiF{;`*EYd2xh{4 zQv#~7%a<>s0%Bq)s)kw7&Xh2FKNd?K#&@i~PBiRZ`9kZotVC~RW!3)UM`Hq?2@3W1 zV6&|~n$>Bh5vTfaHvDj}z~8rHO!>is7Q^#nEjA90(LKNDe2y#VZzJnH>j6h!%*3nT z2p#?LP5NettoOR4#`QQM?|sfP`-zu%y5+NwK(_YxTW6X=wnuFWgq^jdrR~?>DRlLm zYF*~^JTwe3%WUO`slH8148@L#FCzuTb3$tF6pNmZu$0y5c7FRK7EoNn2qkrNXgQQW zEN420h3XAKE;{tOU3bU!yhYu%yE{G=Jv;}I!DGr-@Qc#B=Z}lcF;Hc;BhTNwxil~^ zz|70L^|tKY$B&V_yRLF?6{Oi-EIc&XUMpLaT0XB>XEVPddVVK7zcqEK9*|O*-I(KY zFN=K4i@r2Dlh#+4zbai}SC`->SG%wCzQ#{hk4Cd6%;WG3QoBeAp<1e{(e82|ji@^( zvIy{j-9hrw0~m152~_CheitZH3D?n3&-J4I&c~{}{_x=ymu}f~tje8+P;%Bnvkod` zXJ>Xi7A-jHdQNN(DQz*H$&H;v9B&_%avGLn$fj0NQIXevT&VMdaDQ|a3;MS@r=tHs|I6&_!-(|VgvZBj2|;b=*XZ2w4D zOZ2T<7h~9rA#KuQsz}c>L#{vMtxEa+HOb`4;A;?4XCU0FZL(y2n)@e1CGM3TeEJ-XV^AK_%mTqcJ>eQM4)H{ zY>(HkO5Z8Ridf3#$y$3KDC{W#!LY&8n_t|Oj&nK2!w#?}IX;J2rf#uWGUKY(#S=hT zq1ZdWUm173U5OCHed@{%@DJV49D>KvUVcOt>^8OZT{Xf7z_>K#@w#}qpl;b;rQ!E< z(SU8(QXrgfY4Yo9j<(3M{=3t{T;J`6Y$7U%OC1^!JC!~Guxvzuc%-e6O<#ZfNRp&{ z#X~__nwuphQLEUj$3dJ$iv5q{5sFU2PfR&UBnT#WI#P6_{w3v&dG4^#c6I=SleQx& z>4P^l(4_ICh4>da>YL_TEfb?zbGy2un4bePI@0yk@N$P$tM6#47K&-u49xpQffnO;JwaIFXULiv4uwez7Ffy z*ED>Qztu1w<3{!}?L)=N6XC)(MsK1=j~>y$gaE|Z8k;)P+=vMc)l>8+_%a+g;p3Z7 z<$Of-D@x+*ss}ziUEZXv)t~CGlLL*S6X6}@u%4N2AAfua`SaCUDn21P5Lu{-5rxzD zV)1O#RpIG3AI*vVTVpPTh~}j7Bi&?O)FqE+&*+hhg1tE*9Zw?4fuiwW%^Mr73k6uL=uc%y z*RlgvmclR|{ES%Ih~z?U{#r|?|F??=LIpcxl^^_`{I4%kpWVs3;~^)8Co4UVOW=0Q z;bFo(P2Jt+GvXCL*%f5-Z-npQNGIR^dZE3-lwOWiKl1Rad9aN(M*#6rzOP z#J5Kq6v!iGz2z~i)!RVpsYyycA7?K~wB#?I?7pJY^=XhWRIAu>uS_=c9=$XaA9{Y* zuK}&COJgL6oGritWS1|8#of^Qp|6}s#i35DF%!WEqUM`9;?#YuU^V}`mGBPqhA#1WL#C-noL(!@? z!~4(F=yLR{-DDs?CtHTKYd!+NjcJEdQ-vocsKesZNzxGD~+yEgz0}(t1sf? zsTVun1@hD$N^X4Bi34=^QQZ57e(=Qc{(3kd5G<6}?((zcWlO+xKlC1TTq3~v=CMhF z;-!s%-`o;R8JamQkqn}RmVM~T`%d@c;^Krnx0r$Sy@1;K^GDLiD5dk!h(pgu@d`lm zeRcl+mG@V8akW{4f`XQZ3zGLLb@JR-0&|uD%s9Ka@I7h1P*88X?E{6nlGnTf4-P9B}G9 zcZYqVfL-HwrON_$6HD_E(*|}-1jf7kv*EJL?3^4vlUAbY?N+)68$#1BHLg=eGD(nt z(6F{2|+gS`}edOj~+c{^55waG59E(asXa zxKPGz({{3VIbVAK+D&}Kf1u23{-dtKZGE-n@GkMocH59ioI)pHoo#@_X!+4!w7FTKRenaN3Lo5{-+r z9nVh>2yfliztqywBCn!?=a}Zb@eGf#{i zY?{-G&F9b$@DjF*vESJs0ZT#-EDdB+&@QvR{KUkB&vjYVcCz|1%BVz@3BQ@9dFH7GroM8g|e$n{9rs#aLVT*LS@LlqIm`}=M|kG2|NiZV4))(nICzm zAnKoC6`6OvHy=ZG`Rv@B(|mh$_i3**_sd&OWTj0A=sr1gI@xL>y~@QUqpnV&uCD(2 z-8=jU*}&KktQH$b4|wHe~P8mIUoIQU6t0w0&X$^P2NP@7I3Rk1}+D}-j?e0YuOQL8pkKS7Cu zkQs2sbMuO!+Jf9gtXppHc9hi}Dsr)kENnD2{oU~BE!M>3;;PN$vy?RNz)S6OKN(zZ2Ya_*YMMXswbCIIF5W#?)DXFYMmT@Te zam9Ms7_alZ{QSbg$fJ0D%8bkU{~l1W2a7uA+bZ zAB4JPsSXhh@7)b&`^!~ibZU43x%c z%Ftl)*NIW{b6w)pm!C)vTnB~|ISA&(Ba5rJff(1WUE_VHq@3RgD4P4z_62!FNfrhx zgwFn^|2bWMhDwNH`q!gxkbgI74h$e>jC_w1CdKNM5pu?6t3rqA4TmuKg+khIZfP<8 z9>dO?*d^pVe<7HdE*Q2xQl3L-$xu#u4#9c4A&}4i%!i0hfC%8Gw4EInYC0|IP9sv= zj4?os4V45LH8mb+13o~*91a%v1%ZYrpC2{4A{5#rBV>)u#?#8-`~jUbnVJ>@$vV6% zQIHjC*N-1RdhYhgZTHDB%>cS+9LRdOogs0G?^VMHXN4IeUTUAEMtzTo zAVi5@A1y_=Gsl&m=^P&8zp6ZVP}A_#qL5WXR{ArrP?7Y)KWDz_Q(|=n2|f9b&wkZ? zXyOQ;EabI&H6|toH7*oINkM^$^4xBrj(+n78^)gdxVZJ3a-H9)mFRf zB^BF3VPat=+WEtQ{RJ5icoge(on1g55Xu0)^repKlVf#ujazGdPCNg}Mxk$U!gp1H zVzDPVit>i$ZSQqZaWnJtgZnI{_4OG}PEJs_>kRG#7d1RwJ|fee?00^vm5|@=kHqj6&1t`(Ev^KOol>*P-KF3kzIN71Q}*+_^Y!aj6bexrmiyAkFI~be9Tr=sm;lnppTn{{I<_3Q}^9vzK%s_R(2=6vuJ{;Jwoj&S?bg z>5=IG@TXHn%jAFJ#Hm+72HBtqGMA`cIY_?bO95=w_bfU)Z+Y+Ac77DE`u-Ja@LAXu zn5ayR>-2YgXizDAmMI$7w>PRb1;)BfqQryiPS+ifdxPv{yRix}ylg~&g`37`AiwVR zn21)ZiJY2(a#R?8YAQsqo_ zfrvy3mf2N8Do)F7ejXkko1AWazWR21yQCX7>w9k-^GX*Yq#{`4qyyjf7YT%0w zr}%`kI5YEXFh|SeaBC(rI~x){W@{8v6pK`FyFM*`N~)RW=YHLu*dxRVqT_S&P=Oyc zy+7#({)o1XTP+|72lR0^JWQWXAwyGxK(@C&Mh!j>926XEtiHa!tDKys0cushXb<~- z&VdWS$M+9GzQ(c>SV(qpZw>h-hn$I5&HUNZAtG7FAWf%$A&1vu+yNWR4Kq$4g><<$qt4-sW_QS@c(l4^`jOaWP)ehWQ>w@j`2<^c8{|j(ORP!`Z)W zap`)J$lbVEx}Whw-~d2>Ypa#cp?+$rZe$fIZh0jo=_gOnAwSeZf^d$XoV>oXb+9KH z(QqM79Rr6g^9{k9N5fdC=cm<1W%6CmBbEU1Bui4vJKLu9gG8qqHiVp?Lj?0Rhvp_5SYXoRT&Gad z!LD-Ph0p_|%z<86XM0>22&JD0)neF?V-1t$9z-3@I`aM=8>8U!%ygqFnCq-Yf)4|l zcQK8k*)lH?WNIi%?{7r%1xle%;-}TQOTC^OuQGjKTquNs$o5wuhExcN43uc7=|L?+ zw^d!U*i{>@U+&4*5ZU_r9>)N!QszM7%;qY-N;)PSDxuW8yLo+J24Muo#mU=20(=oa zgZEor^ZMh`Ml-;Dp#ZLGX^tJ>zx_e9qff7G@CQKxCxuV-8(v+`%v%a5HqWd#ya`|- zx$D3?*Jd+#wSa$C(uvYX{L#?es?~iXY}Xy@uczS7Sqo-+u7vSrcr`jXd2^`z+n0PU zQ11N`=pzu1Na07UDsX<(x(tKKh+QcLI|x(tI;y6R2+U`b*VFlXe^?U6|9Rz3lc_nY6H=_L$&&^ZZnTh6oo>-OS?sxEy} zc5|UNIq5NMBOz(n#$*6nOmBY=;YUTB$#a9ildrV$hwQy?Bzh9;tp<*aFu?61sbwSW zH-cl^H8?w&Nk+6hk(77yxAO-^BM@5T+}#D#A3nrDIXOA?THfyF4wM%t>rQGJr!un3 zAAuL7Sdyfa7ky>$fx zIgCIeS!}*-PoJQ-CVu8#)7m@l`yUpQ;FGU98bj|OhNo)CXJ+4;nFS@RY7d*DF_6q# z7ga_5eO*|*98Lk~I49L&;6V(pLL<(dP$P-5 z?o1hA(~=nfw|ZylU-FfZy8~+ALpbb+yd|WsjS=`5ldZmA3M2OM@sXtU-Dh-I(m)sP zlI}NcZ80DI^-B5L{bIv_F+nVsrzbgDwM{t#4{Uf7uiw0Rp>pR48;BJQ6iBJ`ZIKKv z%YM?}IuMrp?}ig)n;R*E6Uzn~sMq$Uv>F8UCp!*4bA~?b0UMKyG@V%*VfRu^p;CEcwj)|AFu+6*faWOG}a58UZZ1DZG z9xmVUJj}-+x)~GEkCmHR$`ACP&%uO;=X|srB9~jl-y2~6McehOSENc$6Du9FLR~Fn z?r5Px=P_R3NbUL%E$88E$c9r_`C=S z0j~ghZ;s-7-uBTgc{sE;@^af&Zh$lIF{wA5?Xx?K6=^CSp4l1~Ha5b@vAnlW-rYs! zE)S>de^GBO;v}U7t$fbb%p;l}=xT)AkR(+v*yzIb2-hKlIKa z>+2xpGd@FnM@zlSZ@-TQSwelOuTKn;Y>0#IXLH`@P~`o4-)~&a&y;!^fsyg@x>~+q znGM}X0%p!$JaHd!K9R42jI@*uL+?6~2k|9q!y*;?r!-uTBf!VU&&BB;S%dU->af#%7ZfQ|*XyH6RwHDx?;5Ks2A3D5V$ z_cbkLBL>UB^^Ed*CCSW_qDq`Umt+)rC--bHjNy9C;OiXBZ9rBm)9Q=UbNs~Mcj~3=O4u! zHM6#(TY)%2T6z%{8>JcysTo+POBjq~O%%ThQG~|_U zTaAZ$ClbI4RWfX&s%j8bHk3>Zh|I z%=I!96K67Q;5J`sBv)@2_>oM(&ZSC=nfLh{pU=d% z2b!;TaGXYp$no;ltDqlsMdhJiM(7wYY2S&^(K z-ySJiY^tX9X*Us~0O~n{D!#imf)tY7+##Xk+CH9SwVo&TKKsj8r4Lr@hecWjadUsU zF=Bz;#rymlX_boM`uz13P2n0GcZe${Ki*tR;J2Yd8G)Y!cF4U44~U`IOVTGHcf$ag z=hso)!GAL}eYlowy>n(x5s<8CvzWWwSH9=6nbTcXW&eE#VAf| zU|;~m53nQ|QP?;Oov zWV>?31crYIQZ_;l*^Yc|88IUxe_Zq7&ZVXEnBX2phM@VgXN5~Dzh|#w`4!Jx{d8V( zadZ8=sNd%E;)lr~#L$eemj7MUMkqg$i75FlBeb(Pypnu#tt#|~*WV2hhPBma{&!+H zPx#G{^_|2l*5-V`txZBg65R4+%AX94ws;&v8+MUrBE zb@F6^x)iV<4~3$jq#_b|DvFL((Vp3w!1z)6ms1n>!_FPtt5vbYAv6;3y>rmiqynKR z^892;0u=x?JPLKT6(WJW6^=4R|KgetP~!OQXRB-f{{5b^G6^Ud>zx!~p@nk?*AAOLXKIq>J{{ge{9|}eg#ySh z-U99b;~1KUhZI;QMBxBR21vrt%y*#1SCrf&zBT*)JbL2S(40)gDXD0ou;`RRTqP3= zW@l3qmP6gC0N6o7LC+cYKR^A)z;mKMQGtJk6wF+^gqBCNH;b{ws@B;!awj*d8^=S* z(#RZw+rWM2VQ$rkaP7ie zb?=?WYi4z?XtHwvA0JKResTxr-@8_EDEEkp1b~4KH3dHqJ?ytF>r|DZluwL-=M_zr zGWbviOS2Ud-uOv2yJWRfDsqc=S~xg?3lVpzU!OK&kN1uBEZR^6TnN^__3z?#SAP~) zFpVQ?vry+5^}kDh;8I?rBto4>0D%O{-qusa=8V}X%0?C zZD-->`tGX8C3zH22B{rivtQkILTru`F*1&O5Th?)c^wPr`c+nwKv0ZJML$Jm68r1s z<*`5!BuGGXzoNy;gUphVd2`T+UhPTvL*;yF1LOt`xv-`VbaI=g)cIhgx%}_Fu4jWK z<{d=lu30X?P{60)c-Z(6zKs*@X;%#sL!RsT``Cekk+1jD^(~nHjZL**~tVpRnbZG4E*CPMbb+PI1fi!f3RCfHzjI_C$(6HF+>hpFdWzDs#oshdFDf7N76YUaARy+ zQ{_#T6%t1O7ajlqEM^D)$58&O>*}k`Ft@af3=4RPZ)oMesrn$mrCo?(VPOf4h={AY zRCVdnB_qgNi1w_cL^SEM^Puv+^&=oHV`5{SAW0y4lmFRH@&fmFZm>gx9-45|>?O#= z?9e{-pCvoR#l)W64Qp2s2J`L@h-%LY#Ngw{kDcQOvcaUw%FW$<+%q-BVzm(ark1G$Xl{Ee{c=AoZ#&2TG`4-lACr_6(cUgkL@z`JxjEsM z1xfAUEMaAb^&>?^0l`m0O*Y0w%Pn@BiyO5``x$JNcu*hv`0LW2mYF;SojkB(gDS-wt83@2xmZP!=K=`_EE@b6&K1J$R> zr>5@wbh8MGl&f*GCtf$e@|wVl0r(zV8n3H9XzG7+#jo)`Xl0EOdk&%?jN}5?K%L{w z$`SDvbbsgg8M4_EQV-1UtpQ3}=Qg?;gM))$M8uu+oFhK`3u98w`JV*s4Hn960`jHcW&=4Or^wfRV zR}8;->lVU|fLR06WwEmXj@b5eAjs&NiL`z_fDsW{4TSf1V8GRZu8_R4vStg|nJ}2| zdU3w)8{=VML`YLfzT+G}OGZ>E3-W;D_N)w0FAYY#Mxdc3h`42sn^jUO#1uNtNFx64 ztc;H_aCqY0yt$iE<)rHl>I~#QD$oGmR%XA=(=NgUQ4-EHl0bhInQYI5Fu{0oz&N1# z=#iQ4NtfVM#ct+kZ?Ed^>M-e2Pcok4+=rjg6JVPE;|>t^u8A8`0S8z%jsHEVUI0J? z;s?kwvf!lP_t?-|o2p}iW=n8ldak6$=U$1|R|8Q`!1k9B!%boout?WAR6cULn3d^Q zI^)2}3xU$54fmF?&Hlj7XAR6rq z={8zPi6eHLoeVz$Lw0}!jEf^*%^>o$K%G+bbjKVO0K)@e>19F!4Pw*=iMbVMiuZZR zu-Hq<{)`9XE~H3M^4S~IqTzdj4T?1%I56Q!uuuqIL*|E{pAzgZ&F}eb-rKGIOuru* z71ahpg2;i?#mcMQ&G!aK#>(uP;HYdvqf!d=5kkuV0{mg>z$2pkQ7p`6IG+VV7SR8N zmBe7kf!m=_KRPjP%0J*dZhgjj;KF5+kRjeF-{6Nkzm+ zv4t)Z?0--}=)jTi3^dsQ46NJ&SnP(eF*;DCp%JyHtQ#?)+Kt|Z=D@*)3=!7fN)YF* zL*AMj*z~-z1=Gk-8`$`SkjnXk*i#B{BryK|{X4hdskgT-6P)2y;bLIwkPZ0m-Me+J zdxUb zR@Kn$AOyRq3R-s906w`@ z-#+t7mbFBO3l@H!JLntlf3$A_k@7?T2N&piAASBPD1q}>x6Br{*7U82KF?d0{)xbEiwa zplNnm`N_V!F+mq4_U8he@gjIwP^4Xd)c$AxH0MkU2Ae)!Q~3NM$vIxo1>XOAZeh~Q z*Be>${V%EHlcWFR+eFMGK5R@58vc4GgG=>L&YqhmEOzUG75Ex2z<6p}G4_02S zeF0EYa|QmVHdRAb-+?%VeExIMV)DPrVXNZuANn9Dp|*x}9gHWdpVA1}(txWh63(2# zMhtznG3xu;=PeR_|HbyV^p%2!hKBSKe)Ncf4<3L?-vj-sp;$J2s_>RYnvny($kM-A#@qZ(qvZW-m89G$N5)wOLCiadyax`X{tfP*9M8{SoOC6Vwp_ zX$aB#`Fsu?A<0ie1eSJoz6tQ9(DXY+LYKK1MEJv%jAU>Vga8qVfjZkSuVWDqpaS=c z=Ujx~`*m-KABZG?n)W|G@aF^17)sDf2hr+LbtS>Jc2l4l2-Y1DgQUo7+=#ur*HtGc zSYWYwOTaQ~!EI-ic(}&syHyCO)I4Iu^ zkDi`cl$|h(mkMinnRk%77s!=EKDf2+@mc;?uNFW&?bV@*323A&$Q=@bDNgDdGRz#6u|24*hszj^`ll{nzFv%mA_4_llf z`@-TP@+!nw2=JkIt$4EfV_uKv&O!Q&&HbtJxBC-Y?`*8CKZB8J1-yIl2L7jenspBF z7XOnurgZ2T*YSd7P1iE?9bBxn^L?A>3qFvaU;xbP87?s7IoQqrdsab3|+Z4|8$0{saxic}FUcE@ze29dWIy}+XO@=mPziT@MQZAPutC6$qRJv8+6Mn6+mt1Lk?hJYkRwB)umiI=CgUnxN&D@bC42J^PG3c69b{t11#cGJd@S4N9Bz}Ze1W?e8*QYd(VYGm8zVaXvA5_Z6E^pc zCxH)={igqBea%kq5kAq`Jv!!G625b%=G)uSi{>P+jw7=OCc< z5g3mve2#9y>QpI5{Y6|t;@8b8F`yj$Rp&1O z_RygjFW4iVuye`b+gANwW&kNlYHMqUFjpKLf3f!en3Y^qj)W5bbHCdH0+MNci1Wm5 z(`|1Fw`3f57XXFVD)d{!f!`e z%mJ=}T@?`x!O}M@VCq?rlEPelw3wJ6?){JG0;G}(n+2d<3xZq05}vDHMMX|d?g{t~ zi6Xy#X@ddscm7P1CKM7}CNz9RoGFXyNV*H)eYt_!KiN_!XwL_hlG zpY2_7K4G{s+rI)MFz?I z({5&IHj-3K@FXFMDMYtsC`9A9LoN+Xh+Xb=BI)dwm9*AW3o$A#RcnC%FGh&cE3;9A zqA`Yk{R(YtlmsK_?s7HsKtV4ZS!rphug2xDwaL4yYo+lyCh-kJ-%QWj25BKUzSI4JO>fP?O1#pYA61S7^jaK8@@Lqpx-+{TLpKB>`0Zu3sEWg_RX zlgD;xRGS(4E%~{U)h@dfkp@>~l#MBGI*i{_;D>(3uN+jH&+-iG7y)`it%Zkz?F&CH zVCtYc-oOiKb$=n7wlaokRc}?t{MqeeXV-TULql%4$}g7}-g78zF%s#<5k!SV-u-pe zW4|*rH=qTpb7+Xl-rinZWNB-r2|5m?O-t#DLumrWM^t?XPR#dFo(onC%GH`` zcqp3RwoyOd%L4n)*{rVwR9t$5i~TC;u@g9=gvI23_qV3v?zf)b|X&9#T$#ILZj zVnI7^Gnhj91_oY?j_Pn4*3v^)eT#k@`;{xKP)VQwAPib0KRZs|-?1>rk8J14ZtWvX z=cS1#R8z>zwK-c5$Mf-@>(;*~(AX92VlGYa+Vt^Tt|MYPT3Q67A+`-v00I+aSm!G_=t5n)HPA)D1#BGJzj8 zc@~W@*`$``hVuz)F{z@@bg!&Q*h`>-yjJT^Y3a^I(?rP5$62Ke(oN4=f{5tlbDT%(O>a9w z$LfWmA}t)69D+4*I^MMsJ3rYu1d9omK5q2LMAgL$DC{U#DW^v!x(S1VtEyAEZAo5S zFE<9}BATd+xcNd^WCsYTlCf5mjEz|UzIA#%h;DVVGO|n({d-b{d1*XfDsU+D-v0GN zQkR+6=*Ps-qsrlX^@!l4W$+XSfzu+s{;NZIZ?6Z*ljn(-klROk&TYwzXIss`M)zSA zY6sxp;$ouyukNndalFC#U->Z^f*_&{54fC&*AZ3 zC4}Uo{KQxLjLwbm0qBj^Fw&dvjA?&lr1}O19lmxkq`E5?U9M()H@+Yfa>$F7Mg&!p z`pMlw73F1NehJR=6Qhy)aCgXpXm@UY9z>u)@Jwt6N{|KQ|I{>0LzBw2+?s@Pue@P` zuRfx>@eV2+Z;?X8M(^pJ&J|l`mb}x}#4>t(?0c~KkMO%|v=RBTRS>ZN@4wRNx6cgF z;ug3FLjWd)N8Q_%$%qZ!kPC;Fq)%Y(bH?Q~tV56t4_A$T;NBH_TW%Vns_XKKV|_lk z_$6F9ZRklW>0)>PS?8PIoo}AG?ufGK%?(ISKD$JbP9{9mN#U@sr!?qeE(Q3ZXfb1^ zD15Zq^KBQ&v#feja||OQcuE})AciP-*R=7_-NUEoufua@F?1)}-``KXV#NgQfSQ`x zXTC?q@@Wu_f$A*{lXbfFI4FT>e4LkQ=*YK?_c|_rsF~B^(g$? zg>sM*oSM#A$Nnu)5$QMQrWv`*X|X zCvaqQ)43J|hZalV4$(DbaeO*CVcRL4&W6rjEOBR8Fnr4WT%sy$Cjt_obWV;(z`dlO zpY`FI$JF^;TWPykHc(>tqvO_fhlt&rHMvKlTm{Pam@HTMwvA?9_r8WzuQSNUR*%6AKB1aTV-5-_K^a@ntaZ_pPrr$fz}+*1YOb@&hgTXRwmh7 zV=H*{^<+QRt4wO$86QR-T}!p<*11N1&h<(KG}5ix^VG(o;sgn z4R`USs6BVuLNoC(sdh>l=35urn+34q5M$z^=T7DiZ9GxG;%|TgYCnDY)HTIOPo1}F zu{jWj4Reb!kTob~^E7MgGe?zncqWzWlhGXSf3Dd)u27@`H4FVL$E8bX2tpu3A#M+7 zb+bSRURhb0*~W4DbA030xXDso|1UXo)GOhDH{VcZm75kNbi!|vNKxr5Dg{>tqT5zl ztabJjaoG4mhfzsnb}te%gXskvYcb!-Vwq-f;+!^k4V~a&sbd^mHRNU1FH=OF#TymL z6&cP39e9i52`6^>k)vx9G28?1YvSO*p!^it*bmi=k8`C5y|FHB-rriP9wt?IvPO~) z7GqDDkL@fTt}4IG>NFW@=IPJ3#Ki>(E>!QlHtqjC`0*B}AXlUjI(S5*UhK;rv(F-9 zo~OF!P3yp38C2hT_~cD1ep{Z}?yqEZT>2Gb>8}})?eg0Y&!?JXSvf;UOCyh#*f9vMTU^bwX3q{Y+(oInw+wc}FzyI$0ud+`aWk=&2(7hdWx%|?v2pXub;>amesuLZF2qR)==9#};i z00R)KdJVo$5N4o+2&}Y6jtClXqEjulf4$Stb|R*{={EVOh1n1|6y>Pk?a$}7Rlvl# zUK>Pn5s>$W26@D#2&TuMKG7_Qw9;s?Rj^^2*(2+_sc&wzj4Ko8UdY%t{{5jTIT^U^ zWW_(L$6Ms;S1<(x1n>s4ps}HumwXd9mY3fpOT4e~r7G>%QH{3y+;aU<9@7J&rfJ4Y z3f?Dvpx$LBzS!=XF`U7H(+w%o;`{ zuQ(P3$5zva+!Qn3{^v`TGASA=?g-krnghzyGR126`t?;$9^$b5_b(M}s_$UQC)Lc@ z;UYL@N*xonldy3&N^)gYCC=##4+{?8rtE*03bhT?vax0%El=b6BX<^BLxMr#+w#(W| zJvw{MzEAE!{r)JEKbm;d-8EdwFpRx#W__LI-=9;lSz2iLs8?<n};`mtOrsyTO5ZK@E6xf0f`+dmNSp&N@&9(4VgqpwZc6{&%vzT8-8Pv0kvPgHEOW%Fap}$~75^sU(`W_MgLeXcxX zSfv@t&)PL6)j%WmlpW36EME&D1#_vr9>k=2f!yKb2m-J8GYy}&XpO5GbuC)K32dyw z|177Kjx%ziV@`>e{qo9~8U69)#`t^t{u-P$yJn9V{7DPo{(D|0rPWAoSz)sCtK>pf z_ehNcUPGVaKs#ER(pXhWr_Ry}`Kk4ok+s&OMZ4vhhI4|EEA#0&B);Fv&pBsS4GpQg zd->(4gIiVGv~wY_6rteC-_e<0p76_L_JMhx1Womi!LThDIFVDpHo6P6eQE3^Z!f1gK45YXL3FqLqtB07iS%eK6fZTQhk~p&sqBfYJ zrP?wXw#~ZM{=wmG^x{3a;gp0+7ILRyFAzHCicbJx0Q4e~L<8)9{m2L#sP#U84LjQz zvTc3r7bgSC-_arN=5S=-xdFAk{Y+6lBh10@ zIn8>AN8-ZLm)Tt=|Htp`U5W({J^G$s z-~RM;L5+9t9F8K%mV3c{938al!e;*~m=Ou5iUdV9WS&p}fbrXDzR;hF{9Rf+oB5C8 zFm-7`)ED*i^t_krs1N}K5X)M*Iss#2#>m=nQmzlOJ`bg5J?!UqV3@3T?_L04ATP^O zB!5d}mi|{Fqf!li^fHG{?@1G$wLz2@Bj@nr5fvl_V)TFFBS0jxbSq2~tk&|IfBe<9 z5NT0kILn!P&)hr%D5XRs?cx3ZNqcY&g~NOba)#jN{|uV9vQuIu;932`;d$ zo0yo$E~IG~7({|i-mE@|hMJD90SMkWOtzkfhX*w+t`;gB zWph$yKrGl?IgO3Cz(RfvFa$s0-HsFI_OCjO2=ejpw5(2dXu#(KN>`Vdm8C|ai=W8nAjT(laW^YYvUJD-ywa-Iiv4uSh$`R#jv+TRRVBj~5K zNrJJNZ{Be7^DEDaw96*djGs97a;yJ@rz61UNV3C-*OEjzi6mg4T4rl&i->F{CQQw7 zxCx**)JYkbv*1Dmsd2qf+^6PdY4FY~F$`K1O03EC+Hv9Ff&-7DA3&sm39 zG2l1xkOq-rp~?^R;CRlRV(6oNFMCr{|Ga+uHhQH{2oYSscn3aI_|OqG`|J63v2DHQ z0g?luzsKNHD-;s77=y>BI9eoM5EA-W+$l*LSO*^NFo@w)Pt2T612f|XAW2arbe$6G zJsOE(eWcem>p-AP0>U2xks*lariJhId4nz=(JNQR8Xhx)IQwN+W+J`Y1QacY6EtD} zEj3<5X5KeUMT78F%P0C+3=fQm$I|-d<1hhAmSFv7gP{ggA__`M)Zn2n1W_Hz|BGdV zgYPOU0PS2b$N?(S4~fE|tCN$HNXi-#8;y8&K^`6ioB)*MWMLDV92N^v4qz+RDKw=- z#4NDyLCsP6p_FJ;6SchTXl`M#z8K0|q(1aGFpCOdg&;dJNy?9{%wbd!a67;^@G&uD zuuqVHHegbq?SBxkF^c}(MoKMF;zf_z5+>>WV^(lL;&0L60hGVacx7aTr;fq@iofWb!o z)+Z;(-O8ZjV%az{=>Lu1ej z)aG{icDA+y)9ph;Ll>=k&chhBR^pxJ3zpj+3YeaefelU!M5F@eO+!NiLgXPSf+)dc z4EOl>;$oyi#EEf0VHKeWAePt!PO18ij>Y$VdNsZhi0uSUJyS?Vv~So zBUXc?GRUz5kFYim;hFJ3Ux(0jK=qHJhg6!iX+C^=1xi^#Rh}qo{P%wVFD2-& literal 0 HcmV?d00001 diff --git a/lectures/cross_section.md b/lectures/cross_section.md index 11754b665..19cd596e6 100644 --- a/lectures/cross_section.md +++ b/lectures/cross_section.md @@ -22,9 +22,7 @@ In addition to what's in Anaconda, this lecture will need the following librarie ```{code-cell} ipython3 :tags: [hide-output] -!pip install quantecon -!pip install --upgrade yfinance -!pip install pandas_datareader +!pip install --upgrade yfinance quantecon pandas_datareader statsmodels interpolation ``` We run the following code to prepare for the lecture: @@ -270,7 +268,7 @@ distribution](https://en.wikipedia.org/wiki/Cauchy_distribution), which is heavy-tailed. (light_heavy_fig1)= -```{figure} /_static/lecture_specific/heavy_tails/light_heavy_fig1.png +```{figure} /_static/lecture_specific/cross_section/light_heavy_fig1.png ``` @@ -461,6 +459,51 @@ plt.show() ``` +## Pareto Tails + +TODO Hi John I added this part with equations you cited below from lecture heavy_tails + +One specific class of heavy-tailed distributions has been found repeatedly in +economic and social phenomena: the class of so-called power laws. + +Specifically, given $\alpha > 0$, a nonnegative random variable $X$ is said to +have a **Pareto tail** with **tail index** $\alpha$ if + +```{math} +:label: plrt + +\lim_{x \to \infty} x^\alpha \, \mathbb P\{X > x\} = c. +``` + +The limit {eq}`plrt` implies the existence of positive constants $b$ and $\bar x$ such that $\mathbb P\{X > x\} \geq b x^{- \alpha}$ whenever $x \geq \bar x$. + +The implication is that $\mathbb P\{X > x\}$ converges to zero no faster than $x^{-\alpha}$. + +In some sources, a random variable obeying {eq}`plrt` is said to have a **power law tail**. + +One example is the [Pareto distribution](https://en.wikipedia.org/wiki/Pareto_distribution). + +If $X$ has the Pareto distribution, then there are positive constants $\bar x$ +and $\alpha$ such that + +```{math} +:label: pareto + +\mathbb P\{X > x\} = +\begin{cases} + \left( \bar x/x \right)^{\alpha} + & \text{ if } x \geq \bar x + \\ + 1 + & \text{ if } x < \bar x +\end{cases} +``` + +It is easy to see that $\mathbb P\{X > x\}$ satisfies {eq}`plrt`. + +Thus, in line with the terminology, Pareto distributed random variables have a Pareto tail. + + ## Exercises ```{exercise} @@ -496,6 +539,11 @@ Use `np.random.seed(11)` to set the seed. ```{exercise} :label: ht_ex4 +(rank_size_fig1)= +```{figure} /_static/lecture_specific/cross_section/rank_size_fig1.png + +``` + Replicate the rank-size plot figure {ref}`presented above `. If you like you can use the function `qe.rank_size` from the `quantecon` library to generate the plots. @@ -533,6 +581,9 @@ Present discounted value of tax revenue will be estimated by 1. multiplying by the tax rate, and 1. summing the results with discounting to obtain present value. +If $X$ has the Pareto distribution, then there are positive constants $\bar x$ +and $\alpha$ such that + The Pareto distribution is assumed to take the form {eq}`pareto` with $\bar x = 1$ and $\alpha = 1.05$. (The value the tail index $\alpha$ is plausible given the data {cite}`gabaix2016power`.) @@ -555,7 +606,6 @@ try to track individual firms given the current distribution. We will discuss firm dynamics in later lectures.) ``` -## Solutions ```{solution-start} ht_ex1 :class: dropdown @@ -596,6 +646,7 @@ plt.show() ```{solution-start} ht_ex2 :class: dropdown +``` Let $X$ have a Pareto tail with tail index $\alpha$ and let $F$ be its cdf. @@ -619,8 +670,6 @@ $$ We know that $\int_{\bar x}^\infty x^{r-\alpha-1} x = \infty$ whenever $r - \alpha - 1 \geq -1$. Since $r \geq \alpha$, we have $\mathbb E X^r = \infty$. -``` - ```{solution-end} ``` From b2840dc28dfd54d3048a2942c2f908e25d9f8367 Mon Sep 17 00:00:00 2001 From: Shu Date: Mon, 24 Apr 2023 14:15:50 +1000 Subject: [PATCH 6/9] fix_label_fig --- lectures/cross_section.md | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/lectures/cross_section.md b/lectures/cross_section.md index 19cd596e6..4988a3220 100644 --- a/lectures/cross_section.md +++ b/lectures/cross_section.md @@ -503,6 +503,10 @@ It is easy to see that $\mathbb P\{X > x\}$ satisfies {eq}`plrt`. Thus, in line with the terminology, Pareto distributed random variables have a Pareto tail. +(rank_size_fig1)= +```{figure} /_static/lecture_specific/cross_section/rank_size_fig1.png + +``` ## Exercises @@ -539,10 +543,6 @@ Use `np.random.seed(11)` to set the seed. ```{exercise} :label: ht_ex4 -(rank_size_fig1)= -```{figure} /_static/lecture_specific/cross_section/rank_size_fig1.png - -``` Replicate the rank-size plot figure {ref}`presented above `. From c4470a9d640909aa1d55aef745584f84f8c96fcc Mon Sep 17 00:00:00 2001 From: Shu Date: Tue, 25 Apr 2023 10:45:42 +1000 Subject: [PATCH 7/9] minor_edit --- lectures/cross_section.md | 174 +++++++++++++++++++------------------- 1 file changed, 85 insertions(+), 89 deletions(-) diff --git a/lectures/cross_section.md b/lectures/cross_section.md index 4988a3220..cca02268f 100644 --- a/lectures/cross_section.md +++ b/lectures/cross_section.md @@ -518,95 +518,6 @@ Replicate {ref}`the figure presented above ` that compares nor Use `np.random.seed(11)` to set the seed. ``` - -```{exercise} -:label: ht_ex2 - -Prove: If $X$ has a Pareto tail with tail index $\alpha$, then -$\mathbb E[X^r] = \infty$ for all $r \geq \alpha$. -``` - - -```{exercise} -:label: ht_ex3 - -Repeat exercise 1, but replace the three distributions (two normal, one -Cauchy) with three Pareto distributions using different choices of -$\alpha$. - -For $\alpha$, try 1.15, 1.5 and 1.75. - -Use `np.random.seed(11)` to set the seed. -``` - - -```{exercise} -:label: ht_ex4 - - -Replicate the rank-size plot figure {ref}`presented above `. - -If you like you can use the function `qe.rank_size` from the `quantecon` library to generate the plots. - -Use `np.random.seed(13)` to set the seed. -``` - -```{exercise} -:label: ht_ex5 - -There is an ongoing argument about whether the firm size distribution should -be modeled as a Pareto distribution or a lognormal distribution (see, e.g., -{cite}`fujiwara2004pareto`, {cite}`kondo2018us` or {cite}`schluter2019size`). - -This sounds esoteric but has real implications for a variety of economic -phenomena. - -To illustrate this fact in a simple way, let us consider an economy with -100,000 firms, an interest rate of `r = 0.05` and a corporate tax rate of -15%. - -Your task is to estimate the present discounted value of projected corporate -tax revenue over the next 10 years. - -Because we are forecasting, we need a model. - -We will suppose that - -1. the number of firms and the firm size distribution (measured in profits) remain fixed and -1. the firm size distribution is either lognormal or Pareto. - -Present discounted value of tax revenue will be estimated by - -1. generating 100,000 draws of firm profit from the firm size distribution, -1. multiplying by the tax rate, and -1. summing the results with discounting to obtain present value. - -If $X$ has the Pareto distribution, then there are positive constants $\bar x$ -and $\alpha$ such that - -The Pareto distribution is assumed to take the form {eq}`pareto` with $\bar x = 1$ and $\alpha = 1.05$. - -(The value the tail index $\alpha$ is plausible given the data {cite}`gabaix2016power`.) - -To make the lognormal option as similar as possible to the Pareto option, choose its parameters such that the mean and median of both distributions are the same. - -Note that, for each distribution, your estimate of tax revenue will be random because it is based on a finite number of draws. - -To take this into account, generate 100 replications (evaluations of tax revenue) for each of the two distributions and compare the two samples by - -* producing a [violin plot](https://en.wikipedia.org/wiki/Violin_plot) visualizing the two samples side-by-side and -* printing the mean and standard deviation of both samples. - -For the seed use `np.random.seed(1234)`. - -What differences do you observe? - -(Note: a better approach to this problem would be to model firm dynamics and -try to track individual firms given the current distribution. We will discuss -firm dynamics in later lectures.) -``` - - ```{solution-start} ht_ex1 :class: dropdown ``` @@ -644,6 +555,13 @@ plt.show() ``` +```{exercise} +:label: ht_ex2 + +Prove: If $X$ has a Pareto tail with tail index $\alpha$, then +$\mathbb E[X^r] = \infty$ for all $r \geq \alpha$. +``` + ```{solution-start} ht_ex2 :class: dropdown ``` @@ -674,6 +592,20 @@ Since $r \geq \alpha$, we have $\mathbb E X^r = \infty$. ```{solution-end} ``` + +```{exercise} +:label: ht_ex3 + +Repeat exercise 1, but replace the three distributions (two normal, one +Cauchy) with three Pareto distributions using different choices of +$\alpha$. + +For $\alpha$, try 1.15, 1.5 and 1.75. + +Use `np.random.seed(11)` to set the seed. +``` + + ```{solution-start} ht_ex3 :class: dropdown ``` @@ -704,6 +636,17 @@ plt.show() ``` +```{exercise} +:label: ht_ex4 + + +Replicate the rank-size plot figure {ref}`presented above `. + +If you like you can use the function `qe.rank_size` from the `quantecon` library to generate the plots. + +Use `np.random.seed(13)` to set the seed. +``` + ```{solution-start} ht_ex4 :class: dropdown ``` @@ -747,7 +690,60 @@ plt.show() ```{solution-end} ``` +```{exercise} +:label: ht_ex5 +There is an ongoing argument about whether the firm size distribution should +be modeled as a Pareto distribution or a lognormal distribution (see, e.g., +{cite}`fujiwara2004pareto`, {cite}`kondo2018us` or {cite}`schluter2019size`). + +This sounds esoteric but has real implications for a variety of economic +phenomena. + +To illustrate this fact in a simple way, let us consider an economy with +100,000 firms, an interest rate of `r = 0.05` and a corporate tax rate of +15%. + +Your task is to estimate the present discounted value of projected corporate +tax revenue over the next 10 years. + +Because we are forecasting, we need a model. + +We will suppose that + +1. the number of firms and the firm size distribution (measured in profits) remain fixed and +1. the firm size distribution is either lognormal or Pareto. + +Present discounted value of tax revenue will be estimated by + +1. generating 100,000 draws of firm profit from the firm size distribution, +1. multiplying by the tax rate, and +1. summing the results with discounting to obtain present value. + +If $X$ has the Pareto distribution, then there are positive constants $\bar x$ +and $\alpha$ such that + +The Pareto distribution is assumed to take the form {eq}`pareto` with $\bar x = 1$ and $\alpha = 1.05$. + +(The value the tail index $\alpha$ is plausible given the data {cite}`gabaix2016power`.) + +To make the lognormal option as similar as possible to the Pareto option, choose its parameters such that the mean and median of both distributions are the same. + +Note that, for each distribution, your estimate of tax revenue will be random because it is based on a finite number of draws. + +To take this into account, generate 100 replications (evaluations of tax revenue) for each of the two distributions and compare the two samples by + +* producing a [violin plot](https://en.wikipedia.org/wiki/Violin_plot) visualizing the two samples side-by-side and +* printing the mean and standard deviation of both samples. + +For the seed use `np.random.seed(1234)`. + +What differences do you observe? + +(Note: a better approach to this problem would be to model firm dynamics and +try to track individual firms given the current distribution. We will discuss +firm dynamics in later lectures.) +``` ```{solution-start} ht_ex5 :class: dropdown From 22fabf1dec8c1884acde89bf84d8e7b110197b5d Mon Sep 17 00:00:00 2001 From: Shu Date: Thu, 27 Apr 2023 14:15:06 +1000 Subject: [PATCH 8/9] minor_edits --- lectures/cross_section.md | 44 +++++++++++++++++++++++++++------------ 1 file changed, 31 insertions(+), 13 deletions(-) diff --git a/lectures/cross_section.md b/lectures/cross_section.md index cca02268f..17caddd61 100644 --- a/lectures/cross_section.md +++ b/lectures/cross_section.md @@ -22,7 +22,7 @@ In addition to what's in Anaconda, this lecture will need the following librarie ```{code-cell} ipython3 :tags: [hide-output] -!pip install --upgrade yfinance quantecon pandas_datareader statsmodels interpolation +!pip install --upgrade yfinance quantecon pandas_datareader interpolation ``` We run the following code to prepare for the lecture: @@ -202,8 +202,8 @@ ax.set_xlabel('returns', fontsize=12) plt.show() ``` -If we look at higher frequency returns data (e.g., tick-by-tick), we often see even more -extreme observations. +If we look at higher frequency returns data (e.g., tick-by-tick), we often see +even more extreme observations. See, for example, {cite}`mandelbrot1963variation` or {cite}`rachev2003handbook`. @@ -220,16 +220,19 @@ frequently. Importantly, there are many examples of heavy-tailed distributions observed in economic and financial settings include -For example, the income and the wealth distributions are heavy-tailed (see, e.g., {cite}`pareto1896cours`, {cite}`benhabib2018skewed`). +For example, the income and the wealth distributions are heavy-tailed +(see, e.g., {cite}`pareto1896cours`, {cite}`benhabib2018skewed`). * You can imagine this: most people have low or modest wealth but some people are extremely rich. -The firm size distribution is also heavy-tailed ({cite}`axtell2001zipf`, {cite}`gabaix2016power`}). +The firm size distribution is also heavy-tailed +({cite}`axtell2001zipf`, {cite}`gabaix2016power`). * You can imagine this too: most firms are small but some firms are enormous. -The distribution of town and city sizes is heavy-tailed ({cite}`rozenfeld2011area`, {cite}`gabaix2016power`). +The distribution of town and city sizes is heavy-tailed +({cite}`rozenfeld2011area`, {cite}`gabaix2016power`). * Most towns and cities are small but some are very large. @@ -304,6 +307,13 @@ def extract_wb(varlist=['NY.GDP.MKTP.CD'], c='all', s=1900, e=2021): return df ``` +```{code-cell} ipython3 +c='all' +s=1900 +e=2021 +wb.download(indicator=['NY.GDP.MKTP.CD'], country=c, start=s, end=e) +``` + ```{code-cell} ipython3 def empirical_ccdf(data, ax, @@ -365,8 +375,13 @@ def empirical_ccdf(data, ### GDP ```{code-cell} ipython3 -df_gdp1 = extract_wb(varlist=['NY.GDP.MKTP.CD']) # gdp for all countries from 1960 to 2022 -df_gdp2 = extract_wb(varlist=['NY.GDP.PCAP.CD']) # gdp per capita for all countries from 1960 to 2022 +# get gdp and gdp per capita for all regions and countries in 2021 +df_gdp1 = extract_wb(varlist=['NY.GDP.MKTP.CD'], s="2021", e="2021")[48:] +df_gdp2 = extract_wb(varlist=['NY.GDP.PCAP.CD'], s="2021", e="2021")[48:] + +# Keep the data for all countries only +df_gdp1 = df_gdp1[48:] +df_gdp2 = df_gdp2[48:] ``` ```{code-cell} ipython3 @@ -458,7 +473,6 @@ fig.tight_layout() plt.show() ``` - ## Pareto Tails TODO Hi John I added this part with equations you cited below from lecture heavy_tails @@ -475,7 +489,8 @@ have a **Pareto tail** with **tail index** $\alpha$ if \lim_{x \to \infty} x^\alpha \, \mathbb P\{X > x\} = c. ``` -The limit {eq}`plrt` implies the existence of positive constants $b$ and $\bar x$ such that $\mathbb P\{X > x\} \geq b x^{- \alpha}$ whenever $x \geq \bar x$. +The limit {eq}`plrt` implies the existence of positive constants $b$ and $\bar x$ +such that $\mathbb P\{X > x\} \geq b x^{- \alpha}$ whenever $x \geq \bar x$. The implication is that $\mathbb P\{X > x\}$ converges to zero no faster than $x^{-\alpha}$. @@ -727,11 +742,14 @@ The Pareto distribution is assumed to take the form {eq}`pareto` with $\bar x = (The value the tail index $\alpha$ is plausible given the data {cite}`gabaix2016power`.) -To make the lognormal option as similar as possible to the Pareto option, choose its parameters such that the mean and median of both distributions are the same. +To make the lognormal option as similar as possible to the Pareto option, choose +its parameters such that the mean and median of both distributions are the same. -Note that, for each distribution, your estimate of tax revenue will be random because it is based on a finite number of draws. +Note that, for each distribution, your estimate of tax revenue will be random +because it is based on a finite number of draws. -To take this into account, generate 100 replications (evaluations of tax revenue) for each of the two distributions and compare the two samples by +To take this into account, generate 100 replications (evaluations of tax revenue) +for each of the two distributions and compare the two samples by * producing a [violin plot](https://en.wikipedia.org/wiki/Violin_plot) visualizing the two samples side-by-side and * printing the mean and standard deviation of both samples. From 96c9fd031a406fe0c466d1fb361aa2ddc9810971 Mon Sep 17 00:00:00 2001 From: Shu Date: Fri, 28 Apr 2023 16:22:00 +1000 Subject: [PATCH 9/9] fix_region --- lectures/cross_section.md | 42 +++++++++++++++++++++++---------------- 1 file changed, 25 insertions(+), 17 deletions(-) diff --git a/lectures/cross_section.md b/lectures/cross_section.md index 17caddd61..fc538afaf 100644 --- a/lectures/cross_section.md +++ b/lectures/cross_section.md @@ -301,19 +301,24 @@ TODO Review exercises below --- are they all too hard for undergrads or should we keep some of them. ```{code-cell} ipython3 -def extract_wb(varlist=['NY.GDP.MKTP.CD'], c='all', s=1900, e=2021): - df = wb.download(indicator=varlist, country=c, start=s, end=e).stack().unstack(0).reset_index() - df = df.drop(['level_1'], axis=1).set_index(['year']).transpose() +def extract_wb(varlist=['NY.GDP.MKTP.CD'], + c='all', + s=1900, + e=2021, + varnames=None): + if c == "all_countries": + # keep countries only (no aggregated regions) + countries = wb.get_countries() + countries_code = countries[countries['region'] != 'Aggregates']['iso3c'].values + + df = wb.download(indicator=varlist, country=countries_code, start=s, end=e).stack().unstack(0).reset_index() + df = df.drop(['level_1'], axis=1).transpose() # set_index(['year']) + if varnames != None: + df.columns = varnames + df = df[1:] return df ``` -```{code-cell} ipython3 -c='all' -s=1900 -e=2021 -wb.download(indicator=['NY.GDP.MKTP.CD'], country=c, start=s, end=e) -``` - ```{code-cell} ipython3 def empirical_ccdf(data, ax, @@ -376,19 +381,22 @@ def empirical_ccdf(data, ```{code-cell} ipython3 # get gdp and gdp per capita for all regions and countries in 2021 -df_gdp1 = extract_wb(varlist=['NY.GDP.MKTP.CD'], s="2021", e="2021")[48:] -df_gdp2 = extract_wb(varlist=['NY.GDP.PCAP.CD'], s="2021", e="2021")[48:] -# Keep the data for all countries only -df_gdp1 = df_gdp1[48:] -df_gdp2 = df_gdp2[48:] +variable_code = ['NY.GDP.MKTP.CD', 'NY.GDP.PCAP.CD'] +variable_names = ['GDP', 'GDP per capita'] + +df_gdp1 = extract_wb(varlist=variable_code, + c="all_countries", + s="2021", + e="2021", + varnames=variable_names) ``` ```{code-cell} ipython3 fig, axes = plt.subplots(1, 2, figsize=(8.8, 3.6)) -empirical_ccdf(np.asarray(df_gdp1['2021'].dropna()), axes[0], add_reg_line=False, label='GDP') -empirical_ccdf(np.asarray(df_gdp2['2021'].dropna()), axes[1], add_reg_line=False, label='GDP per capita') +for name, ax in zip(variable_names, axes): + empirical_ccdf(np.asarray(df_gdp1[name]).astype("float64"), ax, add_reg_line=False, label=name) plt.show() ```