Last active
October 28, 2023 22:04
-
-
Save rueycheng/5c3770df0058f6ee5d87f76a555521ad to your computer and use it in GitHub Desktop.
Python implementation for bootstrapping sampling distributions of Krippendorff's Alpha
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| """ | |
| Bootstrapping sampling distributions of Krippendorff's Alpha | |
| A fast python implementation based on `fast-krippendorff` github repo and | |
| Krippendorff's method. | |
| https://github.com/pln-fing-udelar/fast-krippendorff | |
| http://afhayes.com/public/alphaboot.pdf | |
| """ | |
| import numpy as np | |
| from typing import Any, Callable, Iterable, Optional, Sequence, Union | |
| from krippendorff.krippendorff import ( | |
| _distance_metric, | |
| _distances, | |
| _random_coincidences, | |
| _reliability_data_to_value_counts, | |
| ) | |
| from krippendorff import alpha | |
| def prepare_value_counts_and_domain( | |
| reliability_data: Optional[Iterable[Any]] = None, | |
| value_counts: Optional[np.ndarray] = None, | |
| value_domain: Optional[Sequence[Any]] = None, | |
| level_of_measurement: Union[str, Callable[..., Any]] = 'interval', | |
| ) -> tuple[np.ndarray, Sequence[Any]]: | |
| """ | |
| Perform sanity check and transformation for two input arguments of `bootstrap` | |
| This part of code is a direct derivation from the preprocessing steps in | |
| `alpha` function. It takes four input arguments, `reliability_data`, | |
| `value_counts`, and `value_domain`, and returns only the latter two in | |
| their revised forms. | |
| See the definition of `alpha` for detailed descriptions about these | |
| parameters. | |
| https://github.com/pln-fing-udelar/fast-krippendorff/blob/main/krippendorff/krippendorff.py | |
| Parameters | |
| ---------- | |
| reliability_data : array_like, with shape (M, N) | |
| value_counts : array_like, with shape (N, V) | |
| value_domain : array_like, with shape (V,) | |
| level_of_measurement : string or callable | |
| Returns | |
| ------- | |
| value_counts : array_like, with shape (N, V) | |
| value_domain : array_like, with shape (V,) | |
| """ | |
| if (reliability_data is None) == (value_counts is None): | |
| raise ValueError("Either reliability_data or value_counts must be provided, but not both.") | |
| # Don't know if it's a list or numpy array. If it's the latter, the truth value is ambiguous. So, ask for None. | |
| if value_counts is None: | |
| reliability_data = np.asarray(reliability_data) | |
| kind = reliability_data.dtype.kind | |
| if kind in {"i", "u", "f"}: | |
| # np.isnan only operates on signed integers, unsigned integers, and floats, not strings. | |
| found_value_domain = np.unique(reliability_data[~np.isnan(reliability_data)]) | |
| elif kind in {"U", "S"}: # Unicode or byte string. | |
| # np.asarray will coerce np.nan values to "nan". | |
| found_value_domain = np.unique(reliability_data[reliability_data != "nan"]) | |
| else: | |
| raise ValueError(f"Don't know how to construct value domain for dtype kind {kind}.") | |
| if value_domain is None: | |
| # Check if Unicode or byte string | |
| if kind in {"U", "S"} and level_of_measurement != "nominal": | |
| raise ValueError("When using strings, an ordered value_domain is required " | |
| "for level_of_measurement other than 'nominal'.") | |
| value_domain = found_value_domain | |
| else: | |
| value_domain = np.asarray(value_domain) | |
| # Note: We do not need to test for np.nan in the input data. | |
| # np.nan indicates the absence of a domain value and is always allowed. | |
| assert np.isin(found_value_domain, value_domain).all(), \ | |
| "The reliability data contains out-of-domain values." | |
| value_counts = _reliability_data_to_value_counts(reliability_data, value_domain) | |
| else: | |
| value_counts = np.asarray(value_counts) | |
| if value_domain is None: | |
| value_domain = np.arange(value_counts.shape[1]) | |
| else: | |
| value_domain = np.asarray(value_domain) | |
| assert value_counts.shape[1] == len(value_domain), \ | |
| "The value domain should be equal to the number of columns of value_counts." | |
| assert len(value_domain) > 1, "There has to be more than one value in the domain." | |
| return value_counts, value_domain | |
| def bootstrap( | |
| reliability_data: Optional[Iterable[Any]] = None, | |
| value_counts: Optional[np.ndarray] = None, | |
| value_domain: Optional[Sequence[Any]] = None, | |
| level_of_measurement: Union[str, Callable[..., Any]] = 'interval', | |
| dtype: Any = np.float64, | |
| confidence_level: float = 0.95, | |
| num_iterations: int = 1000, | |
| return_bootstrap_estimates: bool = False, | |
| sampling_method: str = 'krippendorff', | |
| ): | |
| """ | |
| Bootstrap a distribution of Krippendorff's alpha | |
| This algorithm is based on Krippendorff (2006), with the alpha updated | |
| following the coincidence matrix based implementation as in the | |
| `fast-krippendorff` github repo. | |
| https://github.com/pln-fing-udelar/fast-krippendorff | |
| Parameters | |
| ---------- | |
| reliability_data : array_like, with shape (M, N) | |
| Reliability data matrix which has the rate the i coder gave to the j | |
| unit, where M is the number of raters and N is the unit count. Missing | |
| rates are represented with `np.nan`. If it's provided then | |
| `value_counts` must not be provided. | |
| value_counts : array_like, with shape (N, V) | |
| Number of coders that assigned a certain value to a determined unit, | |
| where N is the number of units and V is the value count. If it's | |
| provided then `reliability_data` must not be provided. | |
| value_domain : array_like, with shape (V,) | |
| Possible values the units can take. If the level of measurement is not | |
| nominal, it must be ordered. If `reliability_data` is provided, then | |
| the default value is the ordered list of unique rates that appear. | |
| Else, the default value is `list(range(V))`. | |
| level_of_measurement : string or callable, default: `'interval'` | |
| Steven's level of measurement of the variable. It must be one of | |
| "nominal", "ordinal", "interval", "ratio", or a callable. | |
| dtype: data type, default: `np.float64` | |
| Result and computation data-type. | |
| confidence_level: float, default: `0.95`, | |
| Designated confidence level. | |
| num_iterations: int, default: `1000` | |
| Number of boostrapping iterations. | |
| return_bootstrap_estimates: boolean , default: `False` | |
| When set true, the function will in addition return the bootstrap | |
| estimates (as the second component in the tuple). | |
| sampling_method: str, default: `'krippendorff'` | |
| Pair sampling methods, can be either `'krippendorff'` (Krippendorff | |
| 2006) or `'random'` (pure random sampling). | |
| Returns | |
| ------- | |
| confidence_interval : np.ndarray | |
| Induced confidence interval as an numpy array of shape (2,), in the | |
| form of (lower bound, upper bound). | |
| estimations: np.ndarray | |
| Bootstrap estimates. | |
| References | |
| ---------- | |
| Hayes, Andrew F. & Krippendorff, Klaus (2007). Answering the call for a | |
| standard reliability measure for coding data. Communication Methods and | |
| Measures, 1, 77–89. | |
| Krippendorff, Klaus (2006). Bootstrapping Distributions for Krippendorff’s | |
| Alpha. http://afhayes.com/public/alphaboot.pdf | |
| """ | |
| assert sampling_method in ('krippendorff', 'random') | |
| # Preprocess input arguments | |
| value_counts, value_domain = prepare_value_counts_and_domain( | |
| reliability_data, | |
| value_counts, | |
| value_domain, | |
| level_of_measurement, | |
| ) | |
| distance_metric = _distance_metric(level_of_measurement) | |
| # Compute the unnormalized coincidence matrix and save it for later | |
| V = value_counts.shape[1] | |
| pairable = np.maximum(value_counts.sum(axis=1), 2) | |
| diagonals = value_counts[:, np.newaxis, :] * np.eye(V)[np.newaxis, ...] | |
| unnormalized_coincidences = value_counts[..., np.newaxis] * value_counts[:, np.newaxis, :] - diagonals | |
| # Compute `o`, `e`, and `d` following the original `fast-krippendorff` implementation | |
| o = np.divide(unnormalized_coincidences, (pairable - 1).reshape((-1, 1, 1)), dtype=dtype).sum(axis=0) | |
| n_v = o.sum(axis=0) | |
| e = _random_coincidences(n_v, dtype=dtype) | |
| d = _distances(value_domain, distance_metric, n_v, dtype=dtype) | |
| # Calculate the resampling probabilities based on the unnormalized coincidences | |
| o_raw = unnormalized_coincidences.sum(axis=0) | |
| prob = (o_raw / o_raw.sum()).reshape(-1) | |
| S = np.arange(prob.size) | |
| E = np.eye(prob.size) | |
| # Perform the actual bootstrapping process | |
| est = [] | |
| for i in range(num_iterations): | |
| # Loop through all units | |
| # | |
| # For every unit, calculate the number of needed pairs and resample | |
| # judgments accordingly | |
| units = [] | |
| n_pairs = (pairable * (pairable - 1) / 2).astype(int) | |
| for n in n_pairs: | |
| # Sampling approach 1: Krippendorff's method | |
| if sampling_method == 'krippendorff': | |
| idx = np.array([], dtype=np.int64) | |
| while len(idx) < n: | |
| draw = np.random.choice(S, size=(n, 2), replace=True, p=prob) | |
| idx = np.concatenate([idx, draw[draw[:, 0] != draw[:, 1], 0]]) | |
| idx = idx[:n] | |
| # Sampling approach 2: random sample | |
| else: | |
| idx = np.random.choice(S, size=n, replace=True, p=prob) | |
| assert len(idx) == n | |
| sample = E[idx].sum(axis=0).reshape(o_raw.shape) | |
| units.append(sample) | |
| # From the sample, calculate the alpha value | |
| sampled_coincidences = 2 * np.array(units) | |
| o_new = np.divide(sampled_coincidences, (pairable - 1).reshape((-1, 1, 1)), dtype=dtype).sum(axis=0) | |
| a_new = 1 - (o_new * d).sum() / (e * d).sum() | |
| est.append(a_new) | |
| # Compute sample percentiles | |
| est = np.array(est) | |
| a = 1 - confidence_level | |
| res = np.percentile(est, [100 * a/2, 100 * (1 - a/2)]) | |
| # Include bootstrap estimates as necessary | |
| if return_bootstrap_estimates: | |
| res = (res, est) | |
| return res | |
| EXAMPLE_DATASET = ''' | |
| 1 1 1 2 . 2 | |
| 2 1 1 0 1 . | |
| 3 2 3 3 3 . | |
| 4 . 0 0 . 0 | |
| 5 0 0 0 . 0 | |
| 6 0 0 0 . 0 | |
| 7 1 0 2 . 1 | |
| 8 1 . 2 0 . | |
| 9 2 2 2 . 2 | |
| 10 2 1 1 1 . | |
| 11 . 1 0 0 . | |
| 12 0 0 0 0 . | |
| 13 1 2 2 2 . | |
| 14 3 3 2 2 3 | |
| 15 1 1 1 . 1 | |
| 16 1 1 1 . 1 | |
| 17 2 1 2 . 2 | |
| 18 1 2 3 3 . | |
| 19 1 1 0 1 . | |
| 20 0 0 0 . 0 | |
| 21 0 0 1 1 . | |
| 22 0 0 . 0 0 | |
| 23 2 3 3 3 . | |
| 24 0 0 0 0 . | |
| 25 1 2 . 2 2 | |
| 26 0 1 1 1 . | |
| 27 0 0 0 1 0 | |
| 28 1 2 1 2 . | |
| 29 1 1 2 2 . | |
| 30 1 1 2 . 2 | |
| 31 1 1 0 . 0 | |
| 32 2 1 2 1 . | |
| 33 2 2 . 2 2 | |
| 34 3 2 2 2 . | |
| 35 2 2 2 . 2 | |
| 36 2 2 3 . 2 | |
| 37 2 2 2 . 2 | |
| 38 2 2 . 1 2 | |
| 39 2 2 2 2 . | |
| 40 1 1 1 . 1 | |
| ''' | |
| def get_example_data(): | |
| """ | |
| Get the example reliability data as given in Hayes & Krippendorff (2007, | |
| Table 1) | |
| Returns | |
| ------- | |
| dataset : np.ndarray, with shape (N, M) | |
| Rating matrix of N units and M raters. | |
| """ | |
| rows = [] | |
| for line in EXAMPLE_DATASET.splitlines(): | |
| if len(line) == 0: | |
| continue | |
| row = [int(v) if v.isnumeric() else None for v in line.split()] | |
| rows.append(row[1:]) | |
| return np.array(rows, dtype=float) | |
| def demo(): | |
| """ | |
| Demo function | |
| """ | |
| X = get_example_data() | |
| print('Example Data\n\n', X, '\n') | |
| # Input reliability data is assumed to be of the shape (raters, units) | |
| alpha_value = alpha( | |
| reliability_data=X.T, | |
| level_of_measurement='ordinal', | |
| ) | |
| print(f"Krippendorff's alpha: {alpha_value:.6f}") | |
| ci, est = bootstrap( | |
| reliability_data=X.T, | |
| level_of_measurement='ordinal', | |
| num_iterations=10000, | |
| return_bootstrap_estimates=True, | |
| sampling_method='krippendorff', | |
| ) | |
| mu = np.percentile(est, 50) | |
| print(f"Bootstrap estimate (Krippendorff's method): {mu:.6f} (95%-CI: [{ci[0]:.6f}, {ci[1]:.6f}])") | |
| ci, est = bootstrap( | |
| reliability_data=X.T, | |
| level_of_measurement='ordinal', | |
| num_iterations=10000, | |
| return_bootstrap_estimates=True, | |
| sampling_method='random', | |
| ) | |
| mu = np.percentile(est, 50) | |
| print(f"Bootstrap estimate (random sample): {mu:.6f} (95%-CI: [{ci[0]:.6f}, {ci[1]:.6f}])") | |
| if __name__ == '__main__': | |
| demo() |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| { | |
| "cells": [ | |
| { | |
| "cell_type": "code", | |
| "execution_count": 1, | |
| "id": "f261be01-3e6b-4957-b3a0-57ba8ce0edff", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "import matplotlib.pyplot as plt\n", | |
| "%matplotlib inline\n", | |
| "\n", | |
| "import numpy as np\n", | |
| "\n", | |
| "from krippendorff import alpha\n", | |
| "from bootstrap_alpha import bootstrap, get_example_data" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 2, | |
| "id": "6d97cba9-9d64-4820-a32b-9b67181a2dd9", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "X = get_example_data()\n", | |
| "\n", | |
| "# Calculate the empirical alpha value\n", | |
| "alpha_value = alpha(reliability_data=X.T, level_of_measurement='ordinal')" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 3, | |
| "id": "e703cd9c-93f9-492b-9a78-c10b716dd07d", | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtoAAAEWCAYAAABYLDBhAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAAArLElEQVR4nO3de5hldX3n+/dHEAUEmktLsBtsVNR4iYqtkqgZFI0KUTgninhFJOmYYDKOZiIaLySjc3CeJKhHB4My2nhDvEXG24h4mdEjYHMRFKO00Eg3jbQIqIAi5nv+WL+C3UV19a7LqtpV9X49z372utd3rb3rW99a6/dbK1WFJEmSpNl1j/kOQJIkSVqMLLQlSZKkHlhoS5IkST2w0JYkSZJ6YKEtSZIk9cBCW5IkSeqBhbaGluSkJB9qwwck+WWSHeY7rulK8oEkb5mlbd15bIZY9tAkG6ew7Scn+cH0o+tHkpcl+cYM1t83yf9O8osk/5TO+5PcmOSC2YxV0vDM9ZNua+hcPxeSbEjytPmOQ9tmob0AJHlSkv8vyc1Jfpbkm0keN58xVdWPq+o+VfXb+YxjWDMtCudTVf2fqnrIfMfRgzXAT4Hdq+o1wJOApwMrq+rx8xqZNA/M9TO3kHO9Fqcd5zsATS7J7sBngb8AzgJ2Ap4M/Ho+49LcSLJjVd0x33HMpiQBAtwfuLzuemrW/YENVXXLvAUnzRNz/eK3GPO5ts8z2qPvwQBV9dGq+m1V3VZVX6qqSwGSPDDJV5LckOSnST6cZNnYyu2y0n9OcmmSW5Kc3i7Zf6Fdsv9ykj3bsquSVJI1Sa5NsjnJ30wU1MCyO7bxryX5L+0MzC+SfCnJPgPLvzTJ1S3ONw5zuatdovt4kg+1bV6W5MFJXpfk+iTXJPmjgeX3aPu3OcmmJG9JskOS3wXeA/x+uwR608CP2TPJ59r2z0/ywIHt/UGSb7ezS99O8gcD8w5M8vW23jnAPkxTkr9OcnmSlWnNSpK8Nsl1wPszrqlJO3ava+vc2Jpb3LvNG1v/9e37sCHJiwbWvVeSf0zy4yQ/SfKeJDuPW/c17fhuTnLcwLp7Jzk7yc/TNe14IAO2c7y+luStSb4J3AqcARwL/G37TP4ceN/AZ/T30z2e0gJlrl9kuX4b+XzPJJ9NsqXl788mWTmwzvaO70sGju/fjft590ry9vaZXtuG7zUulr8dyO9HJTk8yQ/TXUF5/bD7puFZaI++HwK/TbI2ybPGEuWAAP8PcD/gd4H9gZPGLfMndJfkHww8G/gC8HpgOd134K/HLf8U4CDgj4DXbi9JDnghcBxwX7qzMX8DkORhwH8HXgTsB+wBrBhym88GPgjsCVwM/K8W8wrgH4B/GVj2A8AdwIOAx7T4/7Sqvg+8AvhWuwS6bGCdY4C/b9tfD7y1xbwX8DngncDewD8Dn0uyd1vvI8CFdEn3v9AVjXdqf+xeuL2dS/Im4GXAf6iqsWL6d4C96M7wrtnGqi8CnkFX7D4YeMPAvN9pca1ocZ2WZKzpyclt+UfTHacVwJvGrTv2+RwPvHvgO/du4Fd0n+HL22tsP7Z3vABe0vZnN7rvyYeB/9Y+k39h68/ozdvYb2mxMtcvzlw/Pp/fA3h/Gz8AuA1417h1Jju+p9Ll0vu1eFcOrPd3wCF0+f1RwOO5+9+Ge3NX3n8v8GLgsXRXT96Y5MBJ9kXTUVW+RvxFl1Q/AGykSy5nA/tuY9mjgIsHxjcALxoY/yRw6sD4XwH/2oZXAQU8dGD+fwNOb8MnAR8at+yObfxrwBsG1vtL4Itt+E3ARwfm7QLcDjxtO/t9EnDOwPizgV8CO7Tx3VoMy4B96S6x7jyw/AuAr7bhlwHfGLf9DwDvGxg/HPi3NvwS4IJxy3+rbeeA9jnsOjDvI2PHZojP81BgE11C/wawx7h5twP3Hjdt47jP9BXj4v7RwLLjYzsLeCPdH+pbgAcOzPt94KqBdW8b+0zbtOvpEvcOwG/GfTf+69gxnex4DXw//mGC4/+WgfG7fUa+fC2lF+b6sfHFlOu3yucTLPNo4MaB8e0d3zMH5u06eHyBHwGHD8x/Bl1zvLFYbpvgmD5hYPkLgaPm+/dgsb1so70AVPdf+ssAkjwU+BDwduAFSfYF3kH33+hudP8t3zhuEz8ZGL5tgvH7jFv+moHhq4FHDhnqdQPDtw5s936D26yqW5PcMOQ2x8f607qrU85t7f0+7WfcE9icZGz5e7D1vkw15qvHLXs13ZmA+9ElxlvGzdt/Oz9r0DK6sxvPr6qbx83bUlW/2s764z+j+w2MTxTb/ejOau0CXDhwjEJXRI+5obZuQzh2TJbT9ekY/3PHTHa8JopZ0jjm+jstply/VT5PsgtwCvBMurPrALsl2WFgf4c9vreMO77j92X834YbJjim2/uOaIZsOrLAVNW/0f13/og26b/S/Vf6yKrane4yUCZee2iDSeQA4NoZbm8zA5e30rUJ3nvbi0/LNXRnOfapqmXttXtVPbzNr0nWnci1dJf2Bh1AdyZ6M117v13HzZuKG4E/pmuz98Rx84aJdbLPaKLYrqW7w8dtwMMHjtEeVTVMYt1Cd2Zn/M8dM9nxGjPVz0Bassz127TQcv34eF4DPITuTPLuwB+26cN8lpsZ+Mxa0T54fMfvy2x8ppohC+0Rl+Sh6TqnrWzj+9NdJjuvLbIb3SW2m5OsAP7zLPzYNybZJcnD6dqJfWyG2/sE8OzW4WQnusuEM/0DsZWq2gx8CfinJLsnuUe6zkP/oS3yE2Bl+/nD+Dzw4CQvTLJjkucDDwM+W1VXA+uAv0+yU5In0V3qnGrMX6Nry/ipJFO9nd0J6TpP7kXXLm/8ZzQW25PpCvqPV9W/07XJOyXJfQGSrEjyjCFi/S3wKeCk9t14GFu3Vdzm8ZrifklLkrl+OAsx14+zG90Jj5ta/p5Kf5RPAH+c7jaQO9G1XR+s4z4KvCHJ8nQdKN9Ed1VE88hCe/T9AngCcH6SW+iS7nfp/iuGrnPHwcDNdB06PjULP/PrdJ1FzgX+saq+NJONVdX36NoHnkn3H/kv6dr+zvZtq15K13Hkcrozxp+g65AD8BXge8B1SX46RMw30BWorwFuAP4W+OOqGlv3hXSfy8/oEuUZg+sn+V4G7vYxyc85h65T4f9McvD2lh/wEbo/NlfStcsbfBjDdXT7fy1dh8NXtLNjAK+l+2zPS/Jz4Mt0Z1eG8Uq6y4rX0Z1pe//AfmzveEmanLl+eAsu1w94O7Az3RXG84AvDrtiO74n0OX/zXT7Pvjws7fQ/WNwKXAZcBFb/23QPEiVV3PVSbIKuAq4Z/V4r88k9wFuAg6qqqv6+jmLVZINdD3svzzBvEPpOuqsHD9PksBcL80lz2hrTiR5drtEuSvwj3T/bW+Y36gkSbPJXC9tzUJbc+VIuqYM19Ldt/WYqqp0D1P45QQvb5wvSQuPuV4aYNMRSZIkqQee0ZYkSZJ6sKAfWLPPPvvUqlWr5jsMSZqyCy+88KdVtXy+45hL5mxJC9l08vaCLrRXrVrFunXr5jsMSZqyJOOfRjfvkvwn4E/pHrJxGd29lfeju13b3nSPaH5JVd2e5F50tzp7LN1t0Z5fVRsm2745W9JCNp28bdMRSRLtISh/DayuqkcAOwDHAG8DTqmqB9Hdt/f4tsrxdI+nfhDdI6XfNvdRS9Jos9CWJI3ZEdg5yY7ALnQPxXgq3QNBANYCR7XhI9s4bf5hSWb1KYCStNBZaEuSqKpNdPc9/jFdgX0zXVORmwYearIRWNGGVwDXtHXvaMvvPX67SdYkWZdk3ZYtW/rdCUkaMRbakiSS7El3lvpA4H7ArsAzZ7rdqjqtqlZX1erly5dU309JstCWJAHwNOCqqtpSVb8BPgU8EVjWmpIArAQ2teFNwP4Abf4edJ0iJUmNhbYkCbomI4e0x2cHOAy4HPgq8Ny2zLHAZ9rw2W2cNv8r5RPQJGkrFtqSJKrqfLpOjRfR3drvHsBpwGuBVydZT9cG+/S2yunA3m36q4ET5zxoSRpxC/o+2pKk2VNVbwbePG7ylcDjJ1j2V8Dz5iIuSVqoPKMtSZIk9cBCW1Pz5UO7l+6uj2Nz6KHdayHwuyGNFn8n7zJMLvV4bZvHZtostCVJkqQeWGhLkiRJPbDQliRJknpgoS1JkiT1wEJbkiRJ6oGFtiRJktQDC21JkiSpBxbakiRJUg8stCVJkqQe7DjfAUgL0aoTP3e3aWc+4AYAjplg3pgNJx/RW0ySpLkz0d+BMeZ6jem10E7yn4A/BQq4DDgO2A84E9gbuBB4SVXdnuRewBnAY4EbgOdX1YY+45PmmolZkqSlo7emI0lWAH8NrK6qRwA7AMcAbwNOqaoHATcCx7dVjgdubNNPactJkiRJC1LfTUd2BHZO8htgF2Az8FTghW3+WuAk4FTgyDYM8AngXUlSVdVzjJIkSXcz2VVIaRi9ndGuqk3APwI/piuwb6ZrKnJTVd3RFtsIrGjDK4Br2rp3tOX3Hr/dJGuSrEuybsuWLX2FL0mSJM1In01H9qQ7S30gcD9gV+CZM91uVZ1WVauravXy5ctnujlJkiSpF33e3u9pwFVVtaWqfgN8CngisCzJWJOVlcCmNrwJ2B+gzd+DrlOkJKlnSR6S5JKB18+TvCrJXknOSXJFe9+zLZ8k70yyPsmlSQ6e732QpFHTZxvtHwOHJNkFuA04DFgHfBV4Lt2dR44FPtOWP7uNf6vN/4rtsyVpblTVD4BHAyTZge7kx6eBE4Fzq+rkJCe28dcCzwIOaq8n0PW1ecLcRy5t33lXbvv2q97xSX3qs432+XSdGi+iu7XfPYDT6BL0q5Osp2uDfXpb5XRg7zb91XTJXJI09w4DflRVV9M1AVzbpq8FjmrDRwJnVOc8uquV+815pJI0wnq960hVvRl487jJVwKPn2DZXwHP6zMeSdJQjgE+2ob3rarNbfg6YN82fGcH9masc/vmgWkkWQOsATjggAP6ileSRpKPYJck3SnJTsBzgI+Pn9ea802pSZ8d2CUtZRbakqRBzwIuqqqftPGfjDUJae/Xt+l3dmBvBju3S5Kw0JYkbe0F3NVsBO7qqA5378D+0nb3kUOAmweamEiS6P/JkJKkBSLJrsDTgT8fmHwycFaS44GrgaPb9M8DhwPrgVuB4+YwVElaECy0pW3w0btaaqrqFsY9kbeqbqC7C8n4ZQs4YY5Ck6QFyUJbkiQteJOdHDlzDuOQBlloS5KkJWusQD/zAdt+qI00XXaGlCRJknpgoS1JkiT1wEJbkiRJ6oGFtiRJktQDC21JkiSpBxbakiRJUg8stCVJkqQeWGhLkiRJPfCBNZIkacGY7AmQo2KyGDecfMQcRqL5ZqGtJW0hJGxJkrQw2XREkiRJ6oGFtiRJktQDC21JkiSpBxbakiRJUg8stCVJkqQeWGhLkiRJPfD2ftKImOhWg2deeQMAh8x1MFqSkiwD3gc8Aijg5cAPgI8Bq4ANwNFVdWOSAO8ADgduBV5WVRfNfdSSNLo8oy1JGvMO4ItV9VDgUcD3gROBc6vqIODcNg7wLOCg9loDnDr34UrSaLPQliSRZA/gD4HTAarq9qq6CTgSWNsWWwsc1YaPBM6oznnAsiT7zWnQkjTiLLQlSQAHAluA9ye5OMn7kuwK7FtVm9sy1wH7tuEVwDUD629s07aSZE2SdUnWbdmypcfwJWn0WGhLkqDrs3MwcGpVPQa4hbuaiQBQVUXXdntoVXVaVa2uqtXLly+ftWAlaSGw0JYkQXdGemNVnd/GP0FXeP9krElIe7++zd8E7D+w/so2TZLUWGhLkqiq64BrkjykTToMuBw4Gzi2TTsW+EwbPht4aTqHADcPNDGRJOHt/SRJd/kr4MNJdgKuBI6jOyFzVpLjgauBo9uyn6e7td96utv7HTf34UrSaLPQliQBUFWXAKsnmHXYBMsWcELfMUmLzUTPTBiz4eQj5jASzQWbjkiSJEk9sNCWJEmSemChLUmSJPXANtrSAmCbPkmSFh7PaEuSJEk9sNCWJEmSemChLUmSJPWg1zbaSZYB7wMeARTwcuAHwMeAVcAG4OiqujFJgHfQPQDhVuBlVXVRn/FJkiSNCvvjLD59n9F+B/DFqnoo8Cjg+8CJwLlVdRBwbhsHeBZwUHutAU7tOTZJkiSpN70V2kn2AP4QOB2gqm6vqpuAI4G1bbG1wFFt+EjgjOqcByxLsl9f8UmSJEl96vOM9oHAFuD9SS5O8r4kuwL7VtXmtsx1wL5teAVwzcD6G9u0rSRZk2RdknVbtmzpMXxJkiRp+vostHcEDgZOrarHALdwVzMRAKqq6NpuD62qTquq1VW1evny5bMWrCRJkjSb+iy0NwIbq+r8Nv4JusL7J2NNQtr79W3+JmD/gfVXtmmSJEnSgtNboV1V1wHXJHlIm3QYcDlwNnBsm3Ys8Jk2fDbw0nQOAW4eaGIiSZIkLSh9P4L9r4APJ9kJuBI4jq64PyvJ8cDVwNFt2c/T3dpvPd3t/Y7rOTZJkiSpN70W2lV1CbB6glmHTbBsASf0GY+WpsnuSypJktSXvs9oS5IkTYknSLRY+Ah2SZIkqQcW2pIkAJJsSHJZkkuSrGvT9kpyTpIr2vuebXqSvDPJ+iSXJjl4fqOXpNFjoS1JGvSUqnp0VY31rzkROLeqDgLO5a7nITwLOKi91gCnznmkkjTiLLQlSZM5EljbhtcCRw1MP6M65wHLxp6RIEnqWGhLksYU8KUkFyZZ06btO/BMg+uAfdvwCuCagXU3tmlbSbImybok67Zs2dJX3JI0krzriBYNe6lLM/akqtqU5L7AOUn+bXBmVVWSmsoGq+o04DSA1atXT2ldSVrohjqjneSRfQciSZod083ZVbWpvV8PfBp4PPCTsSYh7f36tvgmYP+B1Ve2aZKkZtimI/89yQVJ/jLJHr1GJEmaqSnn7CS7JtltbBj4I+C7wNnAsW2xY4HPtOGzgZe2u48cAtw80MREksSQTUeq6slJDgJeDlyY5ALg/VV1Tq/RSZKmbJo5e1/g00mg+9vwkar6YpJvA2clOR64Gji6Lf954HBgPXArcFw/eyNJC9fQbbSr6ookbwDWAe8EHpMuI7++qj7VV4CSpKmbas6uqiuBR00w/QbgsAmmF3DCrAcuSYvIsG20fy/JKcD3gacCz66q323Dp/QYnyRpiszZkjQahj2j/f8C76M7E3Lb2MSquradMZE0Tya728qGk4+Yw0g0QszZkjQChi20jwBuq6rfAiS5B3Dvqrq1qj7YW3SSpOkwZ0vSCBi20P4y8DTgl218F+BLwB/0EZQkaUbM2Rp5PvtAS8Gwt/e7d1WNJWza8C79hCRJmiFztiSNgGEL7VuSHDw2kuSxwG2TLC9Jmj/mbEkaAcM2HXkV8PEk1wIBfgd4fl9BSZJm5FWYsyVp3g37wJpvJ3ko8JA26QdV9Zv+wpIkTZc5W5JGw9APrAEeB6xq6xychKo6o5eoJEkzZc6WpHk2VKGd5IPAA4FLgN+2yQWYtCVpxJizJWk0DHtGezXwsPbIXUnSaDNnS9IIGPauI9+l60wjSRp95mxJGgHDntHeB7g8yQXAr8cmVtVzeolKkjQT5mxJGgHDFton9RmEJGlWnTTfAUiShr+939eT3B84qKq+nGQXYId+Q5MkTYc5W1p8Jntk/YaTj5jDSDQVQ7XRTvJnwCeAf2mTVgD/2lNMkqQZMGdL0mgYtjPkCcATgZ8DVNUVwH37CkqSNCPmbEkaAcMW2r+uqtvHRpLsSHdPVknS6DFnS9IIGLbQ/nqS1wM7J3k68HHgf/YXliRpBqads5PskOTiJJ9t4wcmOT/J+iQfS7JTm36vNr6+zV/V185I0kI1bKF9IrAFuAz4c+DzwBv6CkqSNCMzydn/Efj+wPjbgFOq6kHAjcDxbfrxwI1t+iltOUnSgKEK7ar696p6b1U9r6qe24a9DClJI2i6OTvJSuAI4H1tPMBT6TpWAqwFjmrDR7Zx2vzD2vKSpGao2/sluYoJ2vdV1QNmPSJJ0ozMIGe/HfhbYLc2vjdwU1Xd0cY30t3BhPZ+TdvuHUlubsv/dFwsa4A1AAcccMBUd0WSFrRhH1izemD43sDzgL1mPxxJ0iyYcs5O8sfA9VV1YZJDZyuQqjoNOA1g9erVXgmVtKQM23TkhoHXpqp6O93lRUnSiJlmzn4i8JwkG4Az6ZqMvANY1u5aArAS2NSGNwH7w513NdkDuGFWd0SSFrhhm44cPDB6D7qzJcOeDZckzaHp5Oyqeh3wurb+ocDfVNWLknwceC5d8X0s8Jm2ytlt/Ftt/lfsuyNJWxu2WP6ngeE7gA3A0bMejSRpNsxmzn4tcGaStwAXA6e36acDH0yyHvgZcMw0ty9Ji9ZQhXZVPaXvQCRJs2OmObuqvgZ8rQ1fCTx+gmV+Rdf2W9qmVSd+br5DkObVsE1HXj3Z/Kr659kJR5I0U+ZsSRoNwz6wZjXwF3S3c1oBvAI4mO4WULtNsp5PGZOkuTftnC1Jmj3DttFeCRxcVb8ASHIS8LmqevEQ6449ZWz3Nj72lLEzk7yH7ulipzLwlLEkx7Tlnj/0nmhJ8DKkNJSZ5GxJ0iwZttDeF7h9YPz2Nm1SA08Zeyvw6oGnjL2wLbIWOImu0D6yDUP3lLF3JYm92KXpm+wfkw0ne4fORWxaOVuSNLuGLbTPAC5I8uk2fhR3PXp3Mm9nlp8yJknarunmbEnSLBr2riNvTfIF4Mlt0nFVdfFk6/T1lDEf5ytJk5tOzpYkzb5hO0MC7AL8vKreAWxMcuB2lu/lKWNVdVpVra6q1cuXL59C+JK0pEw1Z0uSZtlQhXaSN9M9tOB1bdI9gQ9Ntk5Vva6qVlbVKroHGXylql4EfJXuKWIw8VPGwKeMSdK0TSdnS5Jm37BntP8v4DnALQBVdS3Tv0XUa+k6Rq6na4M9+JSxvdv0VwMnTnP7krTUzWbOliRN07CdIW+vqkpSAEl2ncoP8SljkjSnZpSzJUmzY9gz2mcl+Re69tV/BnwZeG9/YUmSZsCcLUkjYLtntNu9rz8GPBT4OfAQ4E1VdU7PsUmSpsicLUmjY7uFdrv8+PmqeiRgopakEWbO1lzzib3Stg3bdOSiJI/rNRJJ0mwxZ0vSCBi2M+QTgBe3e2LfAoTuxMnv9RWYJGnazNmSNAImLbSTHFBVPwaeMUfxSJKmyZwtSaNle2e0/xU4uKquTvLJqvqTOYhJkjQ9/4o5W5JGxvbaaGdg+AF9BiJJmrFp5+wk905yQZLvJPlekr9v0w9Mcn6S9Uk+lmSnNv1ebXx9m79q9nZDkhaH7Z3Rrm0MS5JGz0xy9q+Bp1bVL5PcE/hGki/QPan3lKo6M8l7gOOBU9v7jVX1oCTHAG8Dnj/zXZA0VZPd+WXDyUfMYSQab3uF9qOS/JzuLMnObRju6lize6/RSeqNiXlRmnbOrqoCftlG79leBTwVeGGbvhY4ia7QPrINA3wCeFeStO1IkthOoV1VO8xVINIY78kqTc9Mc3aSHYALgQcB7wZ+BNxUVXe0RTYCK9rwCuCa9nPvSHIzsDfw05nEIEmLybD30ZYkLXJV9duqejSwEng83dMlZyTJmiTrkqzbsmXLTDcnSQuKhbYkaStVdRPwVeD3gWVJxq5+rgQ2teFNwP4Abf4ewA0TbOu0qlpdVauXL1/ed+iSNFIstCVJJFmeZFkb3hl4OvB9uoL7uW2xY4HPtOGz2zht/ldsny1JWxv2yZCSpMVtP2Bta6d9D+CsqvpsksuBM5O8BbgYOL0tfzrwwSTrgZ8Bx8xH0JI0yiy0JUlU1aXAYyaYfiVde+3x038FPG8OQpOkBcumI5IkSVIPLLQlSZKkHlhoS5IkST2w0JYkSZJ6YKEtSZIk9cBCW5IkSeqBhbYkSZLUAwttSZIkqQc+sEaSJE1q1Ymfm+8QpAXJM9qSJElSDyy0JUmSpB5YaEuSJEk9sI22JEnSIjVZ+/oNJx8xh5EsTZ7RliRJknrgGW3NC3uwS5Kkxc4z2pIkSVIPLLQlSZKkHth0RNLdbKtpjx1nJEkanme0JUkk2T/JV5NcnuR7Sf5jm75XknOSXNHe92zTk+SdSdYnuTTJwfO7B5I0eiy0JUkAdwCvqaqHAYcAJyR5GHAicG5VHQSc28YBngUc1F5rgFPnPmRJGm0W2pIkqmpzVV3Uhn8BfB9YARwJrG2LrQWOasNHAmdU5zxgWZL95jZqSRptFtqSpK0kWQU8Bjgf2LeqNrdZ1wH7tuEVwDUDq21s08Zva02SdUnWbdmypb+gJWkEWWhLku6U5D7AJ4FXVdXPB+dVVQE1le1V1WlVtbqqVi9fvnwWI5Wk0WehLUkCIMk96YrsD1fVp9rkn4w1CWnv17fpm4D9B1Zf2aZJkpreCm17sEvSwpEkwOnA96vqnwdmnQ0c24aPBT4zMP2lLXcfAtw80MREkkS/99Ee68F+UZLdgAuTnAO8jK4H+8lJTqTrwf5atu7B/gS6HuxP6DE+SdJdngi8BLgsySVt2uuBk4GzkhwPXA0c3eZ9HjgcWA/cChw3p9FKmrFtPTMBfG7CbOmt0G5nNja34V8kGezBfmhbbC3wNbpC+84e7MB5SZYl2c8zJNLomDQpP20OA9Gsq6pvANnG7MMmWL6AE3oNSpIWuDl5MuQMe7BvVWgnWUN3z1YOOOCA/oLWjE1WlEmSJC12vXeGtAe7JEmSlqJez2hP1oO9qjbbg12SpNHgVUhp9vV51xF7sEuSJGnJ6vOMtj3YJUmStGT1edcRe7BLkiRpyfLJkJIkSVIPLLQlSZKkHlhoS5IkST2w0JYkSZJ6MCdPhpQkSdLCMXhf9TMfcAMAx5z4OTacfMR8hbQgWWhrRnzAgSRJ0sRsOiJJkiT1wDPakiQtEV6FlOaWZ7QlSZKkHlhoS5IkST2w0JYkSZJ6YKEtSZIk9cBCW5IEQJL/keT6JN8dmLZXknOSXNHe92zTk+SdSdYnuTTJwfMXuSSNJu86Ikka8wHgXcAZA9NOBM6tqpOTnNjGXws8CziovZ4AnNreNc/G7iwy+JARSfPDM9qSJACq6n8DPxs3+UhgbRteCxw1MP2M6pwHLEuy35wEKkkLhIW2JGky+1bV5jZ8HbBvG14BXDOw3MY2bStJ1iRZl2Tdli1b+o1UkkaMTUe0XYMPOPBSpLR0VVUlqSmucxpwGsDq1auntK4kLXSe0ZYkTeYnY01C2vv1bfomYP+B5Va2aZKkxjPakmbFeVdu+2rHhpOPmOtwNHvOBo4FTm7vnxmY/sokZ9J1grx5oImJJAkLbUlSk+SjwKHAPkk2Am+mK7DPSnI8cDVwdFv888DhwHrgVuC4OQ9Y0pxbNUnTUU+q3J2FtiQJgKp6wTZmHTbBsgWc0G9E2pbJih1Jo8M22pIkSVIPLLQlSZKkHlhoS5IkST2w0JYkSZJ6YGdIAXaskSRJmm0W2pJ65+2gJElLkU1HJEmSpB5YaEuSJEk9sNCWJEmSemAbbUmSRpCd1LXQ2B/n7iy0lxCTtkaRiVmStFhZaC8yFtOSJEmjwTbakiRJUg8stCVJkqQe2HRE0siy/bYWi219l/0ea6lYqvncQluSpHlivxppcRfhFtoLkIlZkiRp9FloS1qQFvMZEEnS4jBShXaSZwLvAHYA3ldVJ89zSPPGs9aSFoKllrfNzdLcWugnVUam0E6yA/Bu4OnARuDbSc6uqsvnNzJJC81CT8wLxWLN2xbT0sKwEHL9yBTawOOB9VV1JUCSM4EjgZFJ2NP9QE3a0uiY7u/jqCTtETMnedscKmmqRqUIT1XN2Q+bTJLnAs+sqj9t4y8BnlBVrxy33BpgTRt9CPCDOQ10avYBfjrfQYwIj0XH43CXpX4s7l9Vy+c7iJkYJm+PUM4e1e+bcU2NcU2NcU3N9uKact4epTPaQ6mq04DT5juOYSRZV1Wr5zuOUeCx6Hgc7uKxWBpGJWeP6vfNuKbGuKbGuKamj7hG6cmQm4D9B8ZXtmmSpNFk3pakSYxSof1t4KAkBybZCTgGOHueY5IkbZt5W5ImMTJNR6rqjiSvBP4X3W2i/kdVfW+ew5qpeb9cOkI8Fh2Pw108FgvcAsvbo/p9M66pMa6pMa6pmfW4RqYzpCRJkrSYjFLTEUmSJGnRsNCWJEmSemChPU1JnpnkB0nWJzlxG8scneTyJN9L8pGB6b9Nckl7LeiOQ9s7DklOGdjXHya5aWDesUmuaK9j5zTwHszwWCyl78QBSb6a5OIklyY5fGDe69p6P0jyjLmNXAvJdHNwkqcM/K5dkuRXSY5q8z6Q5KqBeY+e7bimkxOTPDbJZW2b70ySuYgpyaOTfKsdv0uTPH9gnfk+VhPmzNYx9/y2zY+l66Q7J3GNwHdryrl1mN+jvuJK8vQkF7bv9oVJnjqwztfaNseO133nMK5VSW4b+NnvGVhn6r+LVeVrii+6Tj8/Ah4A7AR8B3jYuGUOAi4G9mzj9x2Y98v53oe5Og7jlv8rus5SAHsBV7b3PdvwnvO9T/NxLJbad4Kus8lftOGHARsGhr8D3As4sG1nh/neJ1+j95ppDh5YZi/gZ8AubfwDwHP7jGvc8kPlROAC4BAgwBeAZ81RTA8GDmrD9wM2A8vm+1i18QlzJnAWcEwbfs9YrpmruObzuzXV3DrVfe0hrscA92vDjwA2DazzNWD1PB2vVcB3t7HdKf8uekZ7eu587HBV3Q6MPXZ40J8B766qGwGq6vo5jnEuDHMcBr0A+GgbfgZwTlX9rB2jc4Bn9hptv2ZyLBaTYY5DAbu34T2Aa9vwkcCZVfXrqroKWN+2J403Wzn4ucAXqurWOYxr0HZzYpL9gN2r6rzq/tKfARw1FzFV1Q+r6oo2fC1wPTBbTzOd9ZzZzi4+FfhEm7SWqR2r2YxrPr5bU82tU93XWY2rqi5u3yuA7wE7J7nXFH/+rMe1LdP9XbTQnp4VwDUD4xvbtEEPBh6c5JtJzksyWETeO8m6Nv2onmPt0zDHAYAk96f7T/orU113gZjJsYCl9Z04CXhxko3A5+nOCA27rgQzz8FjjuHuRdJb22XkU6bxR7+PnLiiDW93mz3ENDjv8XRnBn80MHm+jhVMnDP3Bm6qqju2t80e4xozH9+tk5habp2NnDuTuAb9CXBRVf16YNr7W9ONNw7VRGN24zqwNSn5epInD2xzyr+LFtr92ZHu0uWhdP/xvjfJsjbv/tU94vOFwNuTPHBeIpxbxwCfqKrfzncgI2CiY7GUvhMvAD5QVSuBw4EPJjEXabZNloPHzk49ku4e4GNeBzwUeBzdpf/X9hjfKObECWNqx+qDwHFV9e9t8nwfq1HImZMdr/n4bo1qbp00riQPB94G/PnAOi+qqkcCT26vl8xhXJuBA6rqMcCrgY8k2X2S7UxqFD6AhWiYxw5vBM6uqt+0SzU/pEv6VNWm9n4lXTukx/QdcE+m8vjl8f/dL7ZHN8/kWCy178TxdO0oqapvAfcG9hlyXQlmmIObo4FPV9VvxiZU1ebq/Bp4P1NvutRHTtzUhofZ5mzHRCswPgf8XVWdNzZ9no/VtnLmDcCyJGMP45tODplRXM18fbemmltnI+fOJC6SrAQ+Dby0qu68WjLw+f4C+AhzeLxaE5sb2vQL6a7iPJjp/i5urxG3rwkbw+9I11HlQO5qZP/wccs8E1jbhvehu4SxN10nl3sNTL+CKXY+GJXXMMehLfdQYAPtAUlt2l7AVe147NmG95rvfZqnY7GkvhN0HUhe1oZ/l65dXICHs3WHnSuxM6SvCV4zycED888DnjJunf3ae4C3AyfPdlxtuSnlRO7eAevwOYppJ+Bc4FUTLD+fx2qbORP4OFt3hvzLuYprvr9bTDG3DruvPca1rC3/f0+wzX3a8D3p2ty/Yg7jWk7720PXmXITM/hdnHKC83XnB3Q43RmSH9H9pw/wD8Bz2nCAfwYuBy4b+MX/gzb+nfZ+/HzvS5/HoY2fNFFSAV5O1yljPd0lyXnfn/k4FkvtO0HXu/ubbX8vAf5oYN2/a+v9gCncWcHX0ntNNwe3eavaH897jNvmV9qy3wU+BNxntuNq41PKicDqFtOPgHcxQXHXR0zAi4HftN/Tsdej5/tYTZYz6QqjC9ox/DitIJ/Dz3DevltMI7dOtM25igt4A3DLuO/XfYFdgQuBS+k6Sb6DaZx0mUFcf9J+7iXARcCzZ/K76CPYJUmSpB7YRluSJEnqgYW2JEmS1AMLbUmSJKkHFtqSJElSDyy0JUmSpB5YaGvJSHJUkkry0Da+Ksl3t7POdpeRJM0+c7YWAwttLSUvAL7R3iVJo82crQXPQltLQpL7AE+ie+TqMRPMf1mSzyT5WpIrkrx5YPYOSd6b5HtJvpRk57bOnyX5dpLvJPlkkl3mZm8kaXEzZ2uxsNDWUnEk8MWq+iFwQ5LHTrDM4+meCPV7wPOSrG7TDwLeXVUPB25qywB8qqoeV1WPAr5P9wdBkjRz5mwtChbaWipeAJzZhs9k4kuR51TVDVV1G/ApurMpAFdV1SVt+EK6R+wCPCLJ/0lyGfAi4OF9BC5JS5A5W4vCjvMdgNS3JHsBTwUemaSAHYAC3j1u0drG+K8Hpv0W2LkNfwA4qqq+k+RlwKGzF7UkLU3mbC0mntHWUvBc4INVdf+qWlVV+wNXAfuPW+7pSfZq7fmOAr65ne3uBmxOck+6syOSpJkzZ2vRsNDWUvAC4NPjpn0SeN24aRe06ZcCn6yqddvZ7huB8+mS+7/NQpySJHO2FpFUjb/yIi097TLi6qp65XzHIkmanDlbC4VntCVJkqQeeEZbkiRJ6oFntCVJkqQeWGhLkiRJPbDQliRJknpgoS1JkiT1wEJbkiRJ6sH/D99vYyvwOUbyAAAAAElFTkSuQmCC\n", | |
| "text/plain": [ | |
| "<Figure size 864x288 with 2 Axes>" | |
| ] | |
| }, | |
| "metadata": { | |
| "needs_background": "light" | |
| }, | |
| "output_type": "display_data" | |
| } | |
| ], | |
| "source": [ | |
| "fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(12, 4))\n", | |
| "\n", | |
| "for i, sampling_method in enumerate(['krippendorff', 'random']):\n", | |
| "\n", | |
| " ci, est = bootstrap(\n", | |
| " reliability_data=X.T,\n", | |
| " level_of_measurement='ordinal',\n", | |
| " num_iterations=10000,\n", | |
| " return_bootstrap_estimates=True,\n", | |
| " sampling_method=sampling_method,\n", | |
| " )\n", | |
| "\n", | |
| " mu = np.mean(est)\n", | |
| " ax[i].hist(est, bins=40)\n", | |
| " ylim = ax[i].get_ylim()\n", | |
| "\n", | |
| " ax[i].vlines(x=alpha_value, ymin=ylim[0], ymax=ylim[1], color='red')\n", | |
| " ax[i].vlines(x=[ci[0], ci[1], mu], ymin=ylim[0], ymax=ylim[1], color='orange')\n", | |
| " ax[i].set_xlabel('Alpha')\n", | |
| " ax[i].set_ylabel('Frequency')\n", | |
| " ax[i].set_title(f\"Sampling_method: {sampling_method}\")" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": null, | |
| "id": "2672b42f-db35-49c1-968b-910d834bc400", | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [] | |
| } | |
| ], | |
| "metadata": { | |
| "kernelspec": { | |
| "display_name": "Python 3 (ipykernel)", | |
| "language": "python", | |
| "name": "python3" | |
| }, | |
| "language_info": { | |
| "codemirror_mode": { | |
| "name": "ipython", | |
| "version": 3 | |
| }, | |
| "file_extension": ".py", | |
| "mimetype": "text/x-python", | |
| "name": "python", | |
| "nbconvert_exporter": "python", | |
| "pygments_lexer": "ipython3", | |
| "version": "3.10.5" | |
| } | |
| }, | |
| "nbformat": 4, | |
| "nbformat_minor": 5 | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment