4 "cell_type": "markdown",
7 "# Clustering data using Dhist"
11 "cell_type": "markdown",
14 "Use diagonally cut histogram to model the variability.\n",
16 "Dhist creates a non-uniform histogram which is sampled later to estimate the intercept."
32 "coeff_col_name=\"\"\n",
35 "xbr_file=\"/tmp/dhist_xbr.csv\"\n",
36 "height_file=\"/tmp/dhist_height.csv\""
40 "cell_type": "markdown",
43 "We use a special method to create non-uniform histograms to represent the noise in IO operations.\n",
45 "Unable to install the library properly, I copied the important methods here.\n",
47 "Copied from: https://rdrr.io/github/dlebauer/pecan-priors/src/R/plots.R"
56 "set.seed(seed) #forcing seed for reproductibility\n",
58 "#' Variable-width (dagonally cut) histogram\n",
61 "#' When constructing a histogram, it is common to make all bars the same width.\n",
62 "#' One could also choose to make them all have the same area.\n",
63 "#' These two options have complementary strengths and weaknesses; the equal-width histogram oversmooths in regions of high density, and is poor at identifying sharp peaks; the equal-area histogram oversmooths in regions of low density, and so does not identify outliers.\n",
64 "#' We describe a compromise approach which avoids both of these defects. We regard the histogram as an exploratory device, rather than as an estimate of a density.\n",
65 "#' @title Diagonally Cut Histogram\n",
66 "#' @param x is a numeric vector (the data)\n",
67 "#' @param a is the scaling factor, default is 5 * IQR\n",
68 "#' @param nbins is the number of bins, default is assigned by the Stuges method\n",
69 "#' @param rx is the range used for the left of the left-most bin to the right of the right-most bin\n",
70 "#' @param eps used to set artificial bound on min width / max height of bins as described in Denby and Mallows (2009) on page 24.\n",
71 "#' @param xlab is label for the x axis\n",
72 "#' @param plot = TRUE produces the plot, FALSE returns the heights, breaks and counts\n",
73 "#' @param lab.spikes = TRUE labels the \\% of data in the spikes\n",
74 "#' @return list with two elements, heights of length n and breaks of length n+1 indicating the heights and break points of the histogram bars.\n",
75 "#' @author Lorraine Denby, Colin Mallows\n",
76 "#' @references Lorraine Denby, Colin Mallows. Journal of Computational and Graphical Statistics. March 1, 2009, 18(1): 21-31. doi:10.1198/jcgs.2009.0002.\n",
77 "dhist<-function(x, a=5*iqr(x),\n",
78 " nbins=nclass.Sturges(x), rx = range(x,na.rm = TRUE),\n",
79 " eps=.15, xlab = \"x\", plot = TRUE,lab.spikes = TRUE)\n",
82 " if(is.character(nbins))\n",
83 " nbins <- switch(casefold(nbins),\n",
84 " sturges = nclass.Sturges(x),\n",
85 " fd = nclass.FD(x),\n",
86 " scott = nclass.scott(x),\n",
87 " stop(\"Nclass method not recognized\"))\n",
88 " else if(is.function(nbins))\n",
89 " nbins <- nbins(x)\n",
91 " x <- sort(x[!is.na(x)])\n",
93 " a <- diff(range(x))/100000000\n",
94 " if(a != 0 & a != Inf) {\n",
96 " h <- (rx[2] + a - rx[1])/nbins\n",
97 " ybr <- rx[1] + h * (0:nbins)\n",
98 " yupper <- x + (a * (1:n))/n\n",
99 " # upper and lower corners in the ecdf\n",
100 " ylower <- yupper - a/n\n",
102 " cmtx <- cbind(cut(yupper, breaks = ybr), cut(yupper, breaks = \n",
103 " ybr, left.include = TRUE), cut(ylower, breaks = ybr),\n",
104 " cut(ylower, breaks = ybr, left.include = TRUE))\n",
105 " cmtx[1, 3] <- cmtx[1, 4] <- 1\n",
106 " # to replace NAs when default r is used\n",
107 " cmtx[n, 1] <- cmtx[n, 2] <- nbins\n",
109 " #checksum <- apply(cmtx, 1, sum) %% 4\n",
110 " checksum <- (cmtx[, 1] + cmtx[, 2] + cmtx[, 3] + cmtx[, 4]) %%\n",
112 " # will be 2 for obs. that straddle two bins\n",
113 " straddlers <- (1:n)[checksum == 2]\n",
114 " # to allow for zero counts\n",
115 " if(length(straddlers) > 0) {\n",
116 " counts <- table(c(1:nbins, cmtx[ - straddlers, 1]))\n",
118 " counts <- table(c(1:nbins, cmtx[, 1]))\n",
120 " counts <- counts - 1\n",
122 " if(length(straddlers) > 0) {\n",
123 " for(i in straddlers) {\n",
124 " binno <- cmtx[i, 1]\n",
125 " theta <- ((yupper[i] - ybr[binno]) * n)/a\n",
126 " counts[binno - 1] <- counts[binno - 1] + (\n",
128 " counts[binno] <- counts[binno] + theta\n",
132 " xbr[-1] <- ybr[-1] - (a * cumsum(counts))/n\n",
133 " spike<-eps*diff(rx)/nbins\n",
134 " flag.vec<-c(diff(xbr)<spike,F)\n",
135 " if ( sum(abs(diff(xbr))<=spike) >1) {\n",
137 " counts.new<-counts\n",
138 " diff.xbr<-abs(diff(xbr))\n",
139 " amt.spike<-diff.xbr[length(diff.xbr)]\n",
140 " for (i in rev(2:length(diff.xbr))) {\n",
141 " if (diff.xbr[i-1]<=spike&diff.xbr[i]<=spike&\n",
142 " !is.na(diff.xbr[i])) {\n",
143 " amt.spike<-amt.spike+diff.xbr[i-1]\n",
144 " counts.new[i-1]<-counts.new[i-1]+counts.new[i]\n",
146 " counts.new[i]<-NA\n",
147 " flag.vec[i-1]<-T\n",
149 " else amt.spike<-diff.xbr[i-1]\n",
151 " flag.vec<-flag.vec[!is.na(xbr.new)]\n",
152 " flag.vec<-flag.vec[-length(flag.vec)]\n",
153 " counts<-counts.new[!is.na(counts.new)]\n",
154 " xbr<-xbr.new[!is.na(xbr.new)]\n",
157 " else flag.vec<-flag.vec[-length(flag.vec)]\n",
158 " widths <- abs(diff(xbr))\n",
159 " ## N.B. argument \"widths\" in barplot must be xbr\n",
160 " heights <- counts/widths\n",
162 " bin.size <- length(x)/nbins\n",
163 " cut.pt <- unique(c(min(x) - abs(min(x))/1000,\n",
164 " approx(seq(length(x)), x, (1:(nbins - 1)) * bin.size, rule = 2)$y, max(x)))\n",
165 " aa <- hist(x, breaks = cut.pt, plot = FALSE, probability = TRUE)\n",
167 " heights <- aa$counts\n",
168 " xbr <- aa$breaks\n",
171 " q75<-quantile(heights,.75)\n",
172 " if (sum(flag.vec)!=0) {\n",
173 " amt<-max(heights[!flag.vec])\n",
174 " ylim.height<-amt*amt.height\n",
175 " ind.h<-flag.vec&heights> ylim.height\n",
176 " flag.vec[heights<ylim.height*(amt.height-1)/amt.height]<-F\n",
177 " heights[ind.h] <- ylim.height\n",
180 " end.y<-(-10000)\n",
182 " barplot(heights, abs(diff(xbr)), space = 0, density = -1, xlab = \n",
183 " xlab, plot = TRUE, xaxt = \"n\",yaxt='n')\n",
184 " at <- pretty(xbr)\n",
185 " axis(1, at = at - xbr[1], labels = as.character(at))\n",
186 " if (lab.spikes) {\n",
187 " if (sum(flag.vec)>=1) {\n",
188 " usr<-par('usr')\n",
189 " for ( i in seq(length(xbr)-1)) {\n",
190 " if (!flag.vec[i]) {\n",
192 " if (xbr[i]-xbr[1]<end.y) amt.txt<-1\n",
195 " amt.txt<-amt.txt+1\n",
196 " end.y<-xbr[i]-xbr[1]+3*par('cxy')[1]\n",
198 " if (flag.vec[i]) {\n",
199 " txt<-paste(' ',format(round(counts[i]/\n",
200 " sum(counts)*100)),'%',sep='')\n",
201 " par(xpd = TRUE)\n",
202 " text(xbr[i+1]-xbr[1],ylim.height-par('cxy')[2]*(amt.txt-1),txt, adj=0)\n",
205 " else print('no spikes or more than one spike')\n",
207 " invisible(list(heights = heights, xbr = xbr))\n",
210 " return(list(heights = heights, xbr = xbr,counts=counts))\n",
214 "#' Calculate interquartile range\n",
216 "#' Calculates the 25th and 75th quantiles given a vector x; used in function \\link{dhist}.\n",
217 "#' @title Interquartile range\n",
218 "#' @param x vector\n",
219 "#' @return numeric vector of length 2, with the 25th and 75th quantiles of input vector x. \n",
220 "iqr<-function(x){\n",
221 " return(diff(quantile(x, c(0.25, 0.75), na.rm = TRUE)))\n",
226 "cell_type": "markdown",
229 "## Reading dataframe"
234 "execution_count": 6,
239 "output_type": "stream",
241 " X op msg_size start \n",
242 " Min. : 1130 Length:3500 Min. :1 Min. : 0.2003 \n",
243 " 1st Qu.: 82565 Class :character 1st Qu.:2 1st Qu.:14.3932 \n",
244 " Median :179094 Mode :character Median :4 Median :30.8789 \n",
245 " Mean :174237 Mean :4 Mean :30.3998 \n",
246 " 3rd Qu.:262094 3rd Qu.:6 3rd Qu.:45.7296 \n",
247 " Max. :358579 Max. :7 Max. :63.0828 \n",
248 " duration experiment type index \n",
249 " Min. :1.510e-07 Length:3500 Length:3500 Min. : 2260 \n",
250 " 1st Qu.:1.560e-07 Class :character Class :character 1st Qu.:165130 \n",
251 " Median :1.580e-07 Mode :character Mode :character Median :358189 \n",
252 " Mean :2.322e-07 Mean :348475 \n",
253 " 3rd Qu.:2.270e-07 3rd Qu.:524188 \n",
254 " Max. :1.452e-05 Max. :717158 \n",
257 " Class :character \n",
258 " Mode :character \n",
267 "<table class=\"dataframe\">\n",
268 "<caption>A data.frame: 6 × 9</caption>\n",
270 "\t<tr><th></th><th scope=col>X</th><th scope=col>op</th><th scope=col>msg_size</th><th scope=col>start</th><th scope=col>duration</th><th scope=col>experiment</th><th scope=col>type</th><th scope=col>index</th><th scope=col>origin</th></tr>\n",
271 "\t<tr><th></th><th scope=col><int></th><th scope=col><chr></th><th scope=col><int></th><th scope=col><dbl></th><th scope=col><dbl></th><th scope=col><chr></th><th scope=col><chr></th><th scope=col><int></th><th scope=col><chr></th></tr>\n",
274 "\t<tr><th scope=row>1</th><td>1130</td><td>MPI_Send</td><td>3</td><td>0.200334</td><td>2.57e-07</td><td>grenoble</td><td>exp/exp_PingPong.csv</td><td>2260</td><td>real</td></tr>\n",
275 "\t<tr><th scope=row>2</th><td>1131</td><td>MPI_Send</td><td>3</td><td>0.200338</td><td>1.85e-07</td><td>grenoble</td><td>exp/exp_PingPong.csv</td><td>2262</td><td>real</td></tr>\n",
276 "\t<tr><th scope=row>3</th><td>1132</td><td>MPI_Send</td><td>3</td><td>0.200342</td><td>1.63e-07</td><td>grenoble</td><td>exp/exp_PingPong.csv</td><td>2264</td><td>real</td></tr>\n",
277 "\t<tr><th scope=row>4</th><td>1133</td><td>MPI_Send</td><td>3</td><td>0.200345</td><td>1.56e-07</td><td>grenoble</td><td>exp/exp_PingPong.csv</td><td>2266</td><td>real</td></tr>\n",
278 "\t<tr><th scope=row>5</th><td>1134</td><td>MPI_Send</td><td>3</td><td>0.200349</td><td>1.54e-07</td><td>grenoble</td><td>exp/exp_PingPong.csv</td><td>2268</td><td>real</td></tr>\n",
279 "\t<tr><th scope=row>6</th><td>1135</td><td>MPI_Send</td><td>3</td><td>0.200352</td><td>1.64e-07</td><td>grenoble</td><td>exp/exp_PingPong.csv</td><td>2270</td><td>real</td></tr>\n",
284 "A data.frame: 6 × 9\n",
285 "\\begin{tabular}{r|lllllllll}\n",
286 " & X & op & msg\\_size & start & duration & experiment & type & index & origin\\\\\n",
287 " & <int> & <chr> & <int> & <dbl> & <dbl> & <chr> & <chr> & <int> & <chr>\\\\\n",
289 "\t1 & 1130 & MPI\\_Send & 3 & 0.200334 & 2.57e-07 & grenoble & exp/exp\\_PingPong.csv & 2260 & real\\\\\n",
290 "\t2 & 1131 & MPI\\_Send & 3 & 0.200338 & 1.85e-07 & grenoble & exp/exp\\_PingPong.csv & 2262 & real\\\\\n",
291 "\t3 & 1132 & MPI\\_Send & 3 & 0.200342 & 1.63e-07 & grenoble & exp/exp\\_PingPong.csv & 2264 & real\\\\\n",
292 "\t4 & 1133 & MPI\\_Send & 3 & 0.200345 & 1.56e-07 & grenoble & exp/exp\\_PingPong.csv & 2266 & real\\\\\n",
293 "\t5 & 1134 & MPI\\_Send & 3 & 0.200349 & 1.54e-07 & grenoble & exp/exp\\_PingPong.csv & 2268 & real\\\\\n",
294 "\t6 & 1135 & MPI\\_Send & 3 & 0.200352 & 1.64e-07 & grenoble & exp/exp\\_PingPong.csv & 2270 & real\\\\\n",
299 "A data.frame: 6 × 9\n",
301 "| <!--/--> | X <int> | op <chr> | msg_size <int> | start <dbl> | duration <dbl> | experiment <chr> | type <chr> | index <int> | origin <chr> |\n",
302 "|---|---|---|---|---|---|---|---|---|---|\n",
303 "| 1 | 1130 | MPI_Send | 3 | 0.200334 | 2.57e-07 | grenoble | exp/exp_PingPong.csv | 2260 | real |\n",
304 "| 2 | 1131 | MPI_Send | 3 | 0.200338 | 1.85e-07 | grenoble | exp/exp_PingPong.csv | 2262 | real |\n",
305 "| 3 | 1132 | MPI_Send | 3 | 0.200342 | 1.63e-07 | grenoble | exp/exp_PingPong.csv | 2264 | real |\n",
306 "| 4 | 1133 | MPI_Send | 3 | 0.200345 | 1.56e-07 | grenoble | exp/exp_PingPong.csv | 2266 | real |\n",
307 "| 5 | 1134 | MPI_Send | 3 | 0.200349 | 1.54e-07 | grenoble | exp/exp_PingPong.csv | 2268 | real |\n",
308 "| 6 | 1135 | MPI_Send | 3 | 0.200352 | 1.64e-07 | grenoble | exp/exp_PingPong.csv | 2270 | real |\n",
312 " X op msg_size start duration experiment type \n",
313 "1 1130 MPI_Send 3 0.200334 2.57e-07 grenoble exp/exp_PingPong.csv\n",
314 "2 1131 MPI_Send 3 0.200338 1.85e-07 grenoble exp/exp_PingPong.csv\n",
315 "3 1132 MPI_Send 3 0.200342 1.63e-07 grenoble exp/exp_PingPong.csv\n",
316 "4 1133 MPI_Send 3 0.200345 1.56e-07 grenoble exp/exp_PingPong.csv\n",
317 "5 1134 MPI_Send 3 0.200349 1.54e-07 grenoble exp/exp_PingPong.csv\n",
318 "6 1135 MPI_Send 3 0.200352 1.64e-07 grenoble exp/exp_PingPong.csv\n",
329 "output_type": "display_data"
333 "df = read.csv(text=df_str)\n",
334 "print(summary(df))\n",
339 "cell_type": "markdown",
346 "cell_type": "markdown",
349 "Calculating the intercept for the duration of each message."
354 "execution_count": 7,
358 "df$intercept_estimate = df[[metric]] - df[[coeff_col_name]]*coeff_value"
362 "cell_type": "markdown",
370 "execution_count": 55,
375 "output_type": "stream",
377 "Warning message in hist.default(x, breaks = cut.pt, plot = FALSE, probability = TRUE):\n",
378 "“argument ‘probability’ is not made use of”\n"
387 "4510.135322 5863.669177 7666.843201 2453.222423 1045.499552 347.102832 \n",
388 " 7 8 9 10 11 12 \n",
389 "1448.764591 211.473846 0.000000 0.000000 0.000000 0.000000 \n",
394 " [1] 1.054988e-07 1.210706e-07 1.354500e-07 1.482141e-07 1.828565e-07\n",
395 " [6] 2.538565e-07 4.000706e-07 5.307136e-07 8.678415e-07 1.523917e-06\n",
396 "[11] 2.675975e-06 4.698973e-06 8.251327e-06 1.447928e-05\n",
400 " 1 2 3 4 5 6 7 8 \n",
401 "620.9278 658.0722 690.4405 515.2837 342.9944 157.8873 409.3941 104.0000 \n",
403 " 0.0000 0.0000 0.0000 0.0000 1.0000 \n"
407 "output_type": "display_data"
411 "data = df$intercept_estimate\n",
413 " data = log(data)\n",
416 "dh = dhist(data, plot=FALSE)\n",
422 "cell_type": "markdown",
425 "Saving result to cvs"
430 "execution_count": 56,
434 "write.csv(list(xbr=dh$xbr), xbr_file)\n",
435 "write.csv(dh$heights, height_file)"
439 "cell_type": "markdown",
446 "cell_type": "markdown",
449 "Histogram with vertical lines at the mean of each cluster."
454 "execution_count": 58,
459 "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAC4lBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8RERETExMUFBQVFRUWFhYXFxcY\nGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIkJCQmJiYnJycoKCgpKSkqKiosLCwt\nLS0uLi4vLy8xMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6Ojo7Ozs8PDw9PT0+Pj4/Pz9A\nQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhLS0tMTExNTU1OTk5PT09QUFBRUVFSUlJTU1NU\nVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1eXl5fX19gYGBhYWFiYmJjY2NkZGRlZWVm\nZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29wcHBxcXFycnJzc3N0dHR1dXV2dnZ3d3d4\neHh5eXl6enp7e3t8fHx9fX1+fn5/f3+AgICBgYGCgoKDg4OFhYWGhoaHh4eIiIiJiYmKioqL\ni4uMjIyNjY2Ojo6Pj4+QkJCRkZGSkpKTk5OVlZWWlpaXl5eYmJiZmZmampqbm5ucnJydnZ2e\nnp6fn5+goKChoaGioqKjo6OkpKSlpaWmpqanp6eoqKipqamqqqqrq6usrKytra2urq6vr6+w\nsLCxsbGysrKzs7O0tLS1tbW2tra3t7e4uLi5ubm6urq7u7u8vLy9vb2+vr6/v7/AwMDBwcHC\nwsLDw8PExMTFxcXGxsbHx8fIyMjJycnKysrLy8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS0tLT09PU\n1NTV1dXW1tbX19fY2NjZ2dna2trb29vc3Nzd3d3e3t7f39/g4ODh4eHi4uLj4+Pk5OTl5eXm\n5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u7u7v7+/w8PDx8fHy8vLz8/P09PT19fX29vb39/f4\n+Pj5+fn6+vr7+/v8/Pz9/f3+/v7///9F8m1nAAAACXBIWXMAABJ0AAASdAHeZh94AAAgAElE\nQVR4nO3de5xcZYHm8aMMKmBYHTXLis6oeGNixFEjF8UBxySgIgIaSBBXgdEM4OKKIgyEGVBR\nBy8RFNYri3gZgmZQENFFVhxlFHWQy6AEMVS6w8WI4ZL6f6u6uvvU83S68lRz+nRV9vf9fKhT\np85b5+23c37pk1tTNAE8asVcfwDA9oCQgAoQElABQgIqQEhABQgJqAAhARUgJKACFYX03iW9\nLH7t4iVL9n3hC/Zc+KIlr1zw6vZL+y34m9bjgQv2kYGvWbB3zxO9aK+eh91rxzcLXtbX2/r2\n0gWLu/bGVjsHXrvtIfssOLDqWQd4tY/CtBearfaYqkN62297Hd3UeLDZ/EBRFH+xU/MjxRfb\nL51VfK31eHVxsgz8SfGOntPs+vy+PqrRzqZRvLavt/Xt1cX9XXsPNB6Y3emmMbLtIScWP6h6\n1gcbm6o+ZeSeh2fz7NNdaA837uvefXDJ5FNCqgAh1Y+QCGmWEFJ1CEkQkiGkECEJQjKEFCIk\nQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRB\nSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFI\nhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiG\nkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUGOaQXtELaZYfDCGm2EFJ1CEkQ\nkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCS\nIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIh\npBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGk\nECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQ\nIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAh\nCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGF5iykY/7z4R7+\n0PjTww+Ph/Th4vPtl84sLm09XlX8Dxn44+Ltvc7z8K7P73nYjXY2vy/+tq+39W3/4p6uvU2N\nTbM73TRGtj3khOKaqmf9U+MPVZ8ysnHzbJ59ugttc+Pe7t0HKg/p6Bs39jDSGN24cTyks4rz\n2y99oPh863FNsVIGfq9Y0es8G+ft0fOw29DZ3FIc2Nfb+rZfsa5rb6QxMrvTTaOx7SHHF2ur\nnnV0jla7YXQ2zz7dhTba2NC9e3d5q8OtXQW4tavfdnprR0glQqoDIRHSLCGk6hCSICRDSCFC\nEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIS\nhGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKE\nZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRk\nCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQI\nKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgp\nREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClE\nSIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURI\ngpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiC\nkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQ\nDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAM\nIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwh\nhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGF\nCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUI\nSRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJ\nEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFqg9p\n5Jw3v/E9NzWb93/oyMNPXz91O46QSoRUh2EL6e9PvvV3HzzigeYZJ9925znHPTJlO46QSoRU\nhyEL6b5VrULuXvzrxpJbW1+FDrrBtxPjCKlESHUYspDG/Grp6LWv39J6cvzFvp0YQkglQqrD\nEIZ037Gfba49qv3slNW+bT1ctbDljf/e2JbxkM4sPtXee19xUevxG8U7ZcxVxfKe55i3xzan\n2YqbigNm8rbcvsXtsztBVY4tLp/rD2E4hBfa78qfobcd0h1v+/iW5trl4wHZtvVwwztaDv/F\nPT2MNjbec894SP9YnN9+6QPFF1qPlxcrZeA1xdG9znPPvOf2POwanc2txWv6elvf9ivu7Nob\nbYzO7nTTaGx7yPHF2qpn3ThHq92wcTbPPt2FtrEx0r3b6COkGw6/rPV4XedW7hLfTozi1q7E\nrV0dhu3W7heH/bi9GVlyc7N579IbfTsxjJBKhFSHIQtp8zFfat8LPtA8a+Vt6047YcuU7ThC\nKhFSHYYspBsWj1nT3HTusiNWtS5M344jpBIh1WHIQkoRUomQ6kBIhDRLCKk6hCQIyRBSiJAE\nIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQh\nGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZ\nQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlC\nChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIK\nEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoR\nkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGS\nICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIg\nJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAk\nQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRD\nSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENI\nIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0gh\nQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFC\nEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIS\nhGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKE\nZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSKG5C+k3\nW3rY1Ni8Zct4SOcWX2i/tKr4auvxu8XJMvDfinf0Os+WXZ/f87Ab7WzuLl7b19v69urivq69\nPzb+OLvTTWNk20NOLL5f9aybG5uqPmXknodm8+zTXWgPNbp/pLdsrjyko38+2sNIY2R0dDyk\nVcWn2y+dWnyu9XhZ8S4ZeHWxotd5Ruft0fOw29DZ3Fwc2Nfb+rZvcUfX3thq50Bj20OOK75V\n/bRzs9oNszrtdBfaSGND9+768laHW7sKcGtXv+311o6QJhFSHQiJkGYJIVWHkAQhGUIKEZIg\nJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAk\nQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBCjzKkhb/sbL/yvH4nJqQSIdVh\noEMqrh/bPHT64/qdmJBKhFSHAQ6pKL2o34kJqURIdRjgkG74aLH06La3nnpHvxMTUomQ6jDA\nITWbB/56phMTUomQ6jDQIc0cIZUIqQ4DHdL6I3d7bOcXSf1OTEglQqrDQId0yJ/tf+TYr5KO\n7ndiQioRUh0GOqQnf32mExNSiZDqMNAh7Xz3TCcmpBIh1WGgQ9r7uzOdmJBKhFSHgQ7pxy+5\ndoYTE1KJkOow0CEt2r3Y+Zlj+p2YkEqEVIeBDmnv/Sf0OzEhlQipDgMd0swRUomQ6kBIhDRL\nCKk6j/bPkSbM63diQioRUh0GOqSlY16y0wuO63diQioRUh0GOqRxd+2zpt+JCalESHUYhpCa\n1y/sd2JCKhFSHYYipLt26ndiQioRUh2GIaQtZz6934kJqURIdRjokP5qzAv+vDip34kJqURI\ndRiCkBa86qOb+52YkEqEVIeBDmnmCKlESHUY8JA2rFl9/tr7mn0jpBIh1WGgQ3rkxB3b37Bh\nl7P7npiQSoRUh4EO6ezi4Au+ueZTBxYX9TsxIZUIqQ4DHdLzTuhs38Z3Wu0LIdVvoEN6/JWd\n7eX8gWxfCKl+Ax3SLpd1tl9/Yr8TE1KJkOow0CG94pVjf4D0wAH79TsxIZUIqQ4DHdLlj3nG\n28/4h2N2e+y3+52YkEqEVIeBDqn5tee2f/v7hZf3PTEhlQipDoMdUrN554+u//0MJiakEiHV\nYbBDuutjrYe7T1/f98SEVCKkOgx0SP8xv/3/vLy9mH9rvxMTUomQ6jDQIR307B+1N7989uv6\nnZiQSoRUh4EO6Smf6Ww/xXcR6gsh1W+gQ9rp853tF3bud2JCKhFSHQY6pJcfOPbR3bfXon4n\nJqQSIdVhoENa+5hnHXfaqcuf8ti1/U5MSCVCqsNAh9S8YmH7D2T35A9k+0NI9RvskJrNDT/7\nxQz+gSwhde0RUh0GPaQZIqQSIdWBkAhplhBSdQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBS\niJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKI\nkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQ\nBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAE\nIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQh\nGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZ\nQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlC\nChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIK\nEZIgJENIoVkIad2JS9ub+z905OGnr5+6HUdIJUKqw7CFdM2yc8dCOuPk2+4857hHpmzHEVKJ\nkOowbCFdefcP2yE1ltza+ip00A2+nRhGSCVCqsOwhdRsjoV07eu3tB6Pv9i3E4MIqURIdRjS\nkNYe1X56ymrfth6+/8qWQ3+2oZdG67/xkM4sVrdfeX/xv1qP/1K8U8Z9t1je8zzz9uh5eGvT\nttxUHNDX2/q2b/EbmbYxu9NNI5j12OLyuZh2NszutNNeaPpje1f5M3Qa0vLxgGzberj+zS2H\n3bixh5HG6MaN4yGdVZzffukDxedbj2uKlTLwe8WKXufZOG+Pnofdhs7mluLAvt7Wt/2KdV17\nI42R2Z1uGo1tDzm+WFv1rKNztNoNo7N59ukutNHGhu7du/sO6brOrdwlvp0YxK1diVu7Ogzp\nrd3IkpubzXuX3ujbiUGEVCKkOgxbSKONK5Y2WlfGWStvW3faCVumbMcRUomQ6jBsIR29uO0b\nzU3nLjtiVevC9O04QioRUh2GLaQQIZUIqQ6EREizhJCqQ0iCkAwhhQhJEJIhpBAhCUIyhBQi\nJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIk\nQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRB\nSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFI\nhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiG\nkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQ\nQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBC\nhCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKE\nJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQk\nCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQI\nyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJ\nEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQ\nUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBS\niJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKI\nkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQ\nBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBCcxbSip9t6KXR\n+m88pDOL1e1X/qrY+7DDDtu/eKeM+26xvOd55u3R8/A0bioOmMnbcvsWv+nebTRmd7ppBLMe\nW1w+F9POhtn9JE97oem0d5U/Q8/VV6SJkPiKVBm+IlVneG7tCKlyhFQdQhKEZAgpREiCkAwh\nhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGF\nCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUI\nSRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJ\nEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQ\nkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCS\nIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIh\npBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGk\nECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQ\nIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAh\nCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSKEhCOlFxctbBR1GSJUjpOoQkiAkQ0gh\nQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFC\nEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIS\nhGQIKURIgpAMIYUISRCSIaQQIQlCMoQUIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKE\nZAgpREiCkAwhhYYvpOe3dyYHEtKjQEjVISRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKE\nZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQiJEFIhpBChCQIyRBSaEhDGkdIjwohVYeQBCEZQgoR\nkiAkQ0ghQhKEZAgpREiCkAwhheoL6f4PHXn46esndwmpREh12F5COuPk2+4857hHJnYJqURI\nddhOQmosubX1VemgGyb2CalESHXYTkK69vVbWo/HXzyxT0glQqrDdhLS2qPaj6esbj38cEnL\nYT8f7WGkMTI6uu2QRkevLlb0Os/ovD16HnYbOpubiwP7elvf9i3u6NobW+0c2LDtIccV36p8\n2sFd7aMw3YU20pB515c/Q880pOV9h3RyK6RnPGF0VfHp9kunFp9rPV5WvEsGEtKjQEjVqS2k\n6zq3dpdM7Ae3dh9ohfQXOzU/Unyx/dJZxddaj1fz7biqw61ddWq7tRtZcnOzee/SGyf2CalE\nSHXYTkJqnrXytnWnnbBlYpeQSoRUh+0lpE3nLjti1ejkLiGVCKkO20tIhpBKhFQHQiKkWUJI\n1SEkQUiGkEKEJAjJEFKIkAQhGUIKEZIgJENIIUIShGQIKURIgpAMIYUISRCSIaQQIQlCMoQU\nIiRBSIaQQoQkCMkQUoiQBCEZQgoRkiAkQ0ghQhKEZAgpREiCkAwhhQhJEJIhpBAhCUIyhBQi\nJEFIhpBChCQIyRBSiJAEIRlCChGSICRDSCFCEoRkCClESIKQDCGFCEkQkiGkECEJQjKEFCIk\nQUiGkEKEJAjJEFKIkAQhGUIKzV1IH72wh0+fd8GFF755/vz5z9v9wmPm/337paPmv7v1eOr8\nQ2XgqvkH9TrPhc/cs+dh94nO5rz5r+zrbX17xfzVXXvnn3f+7E43jY9ve8gh80+vetYLzvt0\n1aeMfPIzs3n26S60z5z3qe7dCyoP6XuX9nLmitWXXvrh173ikHeecOk/r/xk+6WPrGy9dOln\nVn5QBl608p96nujE9/Y87C7ubL688h/6elvfTlv5v7v2Przi3NmdbhoXb3vIOSs/W/Wsn1xx\nVtWnjFz8ldk8+3QX2kUr3if7/1p1SL2dt/BHdUwzIC5d+I25/hBqdN3CT8z1h1Cj3yx8/zRH\nCKlyhLT9IqQaEdL2i5BqREjbrzkOCdjeERJQAUICKkBIQAUqDen+Dx15+Onrpz7f6qBpBwyP\nv1vcckjn+VaXM/ni5W89+Pih/w2XdScunXjac7Xdn5bh1bXarS7IVltpSGecfNud5xz3yJTn\n4zPfJAemDhg6yy9rNBrjf8Vtq8uZePE7y65f//Vj5uYvpVXmmmXnTl5aPVfb/WkZWt2r3eqC\nbLVVhtRYcmsrl4NukOejZy97w3tuab/0s3d2H+gePKzecP3k062st+vFY66co4+wSlfe/cOJ\nS6v3ars+LcOra7WTC+q12ipDuvb1W1qPx18sz088+77Nn3vz5uZESBMHugcPqQcXf+xdK1at\nG3u+lfWWL25YfOXfveHEX83dR1qRyUur52q7Py3DbHK1kwvqtdoqQ1p7VPvxlNXdz29ZPNps\nbnnTNc2JkCYOdA8eUve85cM33XTaW/7Qfr6V9ZYv3rT4f95x3+o33TNnH2lFJi+tnqvt/rQM\ns8nVTiyo52orDWl559zdz69ZPOaSGw499JAlhx56wuSB7sHD7I+HXNHe+Hq/v3Tp0l9OvHjT\n4tYd7MOHfWcuP84qlCH1Wu3YiPFPyzArb+3aWgvqudoqQ7qu89Xuku7n1y3ufCXcvH79Nceu\nX7+hPNA1eKgdO/YPFX29m26//fY/TbzYWHxza3vc0K918tLqudrOkM6nZZhpSK0F9VxtlSGN\nLGldL/cuvbH7+W8X/0dr7672S51bu4kD3YOH1O3//FCz+cAhV7Wfb2W95YuPLLus9VPJG6+Z\nw4+1EpOXVs/Vdn9ahtnkaicW1HO1lf7291krb1t32glbmlf8S/n8lGvZR5IAAAR2SURBVHff\n/fA339D+vcNOSJMHJrbD677Dz71r3arlf5pmvc3yxUuO+GnjY8vm5l+gV2a0ccXSRuOBba52\n8tMy1LpWO7mgXqutNKRN5y47YlXrF2Rnv698PvpPh77x3TduZdDEdojd+r5D33zG73usd+LF\nRy56y8Hv6fltLYbA0WO/QvjGtlc78WkZat2rnVhQr9XyV4SAChASUAFCAipASEAFCAmoACEB\nFSAkoAKEBFSAkObMX++xPU3z/ztCmjPnrure+2nfPxLbfsfYCJ2m73MgwydyQHys7x+Jbb+j\nihHI8ImcM+17rr1f8ZNXzXvKm9Y3DyyKYmGzefWr5+204ILWwUV7X/b0lzWbV+zzxKcd0v43\nGJMHXvTSK/fa6UnL75l4R2lyyO/e+ozHP+11vxof0Znmmr2esNvZD5682xP3v7U14kt77TRv\n4ZeaU2bFjBHSnGlf4fvvvte3139lhyObv15aXP/L5nd22OeyK95efLDZfNWezz1vTfOKxxzw\n+Que9V/v6jrwsqe8+P80PrfjwePvKJVDXjr//Ku+8MKnbuqMGJvm6a/8tzsOLl59+rrv7fq3\nzeaXi4PXrHlNsWbKrJgxQpozY1f42P9Cb//dms2j2z8SC57d/k5DS+Y90Drw1dazF//lQ83m\n/33cR7sOLCra/6zp6OK3nXeUJofcW7yn9eSWVXd2RnSmuaHZ/H7x8tbuEbs0m6tetbnZvPfP\njpgyK2aMkObM2BW+c/vZkY/tXNLri3c90PLJ4kfN/R/3YLO5oTi2M7TrwKJd2v+E66LimxZS\nOeTBJz/zO51vlFWG1KqneUtxUuvxpGLif9749L2nzIoZI6Q5M3aFP7P9rH05t//7aTHuq2Nf\npJo/L07rDO06sOhZ7RfWFJ+1kLqG/OAviye//gsPdYfUnuY/i7NajycXG5v3vv8Fu+6wQ7Fo\nyqyYMUKaM1sLacUPxzQ6B35RnNoZ2nWgE9LXi4umhDQ5pPnwlSc9v3jxH6cNaZ8d3nvNz36+\n26Ips2LGCGnOTA1ppDhy4uDYgfuKsW9Vc/vdXQcW7dT+33l/qlhrIXUNGfPx4sLpQrq5OKb1\n5KEnLJoyK2aMkOaMhvTWonUv9pL/srG1e9EpD40feOFTWr+g+VXrBq88sKj1q6Nm86DHj3be\nUZoc8uND29+U+5binM6IqSH9sji92f4jpJdOmRUzRkhzRkM6tTj9K82rd9zzon99345HTRxY\n89i//uLq5zz1rq4Di3Z/zie+fVKxbPwdpckhd83b84Jvf/nlu97SGTE1pAd3/2/f+MGJ++03\n76o/2KyYMUKaMxrSHQt2bO1//2/m7ficsx+aONC8/KU7P/XgXze7Dix67o/32flJx9w/8Y7S\n5JB/P/ipO+528E/GR2zl10jXv2znp/33ey/78yfdZLNixghpyCzi76AOJEIaMoQ0mAhpyBDS\nYCKkISMhfauY9Ik5+4jQRkjD7P6fTxru71k7/AgJqAAhARUgJKAChARUgJCAChASUAFCAirw\n/wBOznRqSUaAyQAAAABJRU5ErkJggg==",
470 "output_type": "display_data"
474 "library(ggplot2)\n",
477 " dh$xbr = exp(dh$xbr)\n",
480 "ggplot(df, aes(x=intercept_estimate)) + geom_histogram(bins=100) +\n",
481 " geom_vline(xintercept=dh$xbr) + \n",
486 "cell_type": "markdown",
489 "Same but in log scale"
494 "execution_count": 59,
499 "output_type": "stream",
501 "Warning message:\n",
502 "“Transformation introduced infinite values in continuous y-axis”\n",
503 "Warning message:\n",
504 "“Removed 59 rows containing missing values (geom_bar).”\n"
509 "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAC6FBMVEUAAAABAQECAgIDAwME\nBAQFBQUGBgYHBwcICAgJCQkKCgoLCwsMDAwNDQ0ODg4PDw8RERESEhITExMUFBQVFRUWFhYX\nFxcYGBgZGRkaGhobGxscHBwdHR0eHh4fHx8gICAhISEiIiIkJCQmJiYnJycoKCgpKSkqKior\nKyssLCwtLS0uLi4vLy8wMDAxMTEyMjIzMzM0NDQ1NTU2NjY3Nzc4ODg5OTk6Ojo7Ozs8PDw9\nPT0+Pj5AQEBBQUFCQkJDQ0NERERFRUVGRkZHR0dISEhJSUlKSkpLS0tMTExNTU1OTk5PT09Q\nUFBRUVFSUlJTU1NUVFRVVVVWVlZXV1dYWFhZWVlaWlpbW1tcXFxdXV1eXl5fX19gYGBhYWFi\nYmJjY2NkZGRlZWVmZmZnZ2doaGhpaWlqampra2tsbGxtbW1ubm5vb29wcHBxcXFycnJzc3N0\ndHR1dXV3d3d4eHh5eXl6enp7e3t8fHx9fX1+fn6AgICBgYGCgoKDg4OFhYWGhoaHh4eIiIiJ\niYmKioqLi4uMjIyNjY2Ojo6Pj4+QkJCRkZGSkpKTk5OUlJSVlZWWlpaXl5eYmJiZmZmampqb\nm5ucnJydnZ2enp6fn5+goKChoaGioqKjo6OkpKSlpaWmpqanp6eoqKipqamqqqqrq6usrKyt\nra2urq6vr6+wsLCxsbGzs7O0tLS1tbW2tra3t7e4uLi5ubm6urq7u7u8vLy9vb2+vr6/v7/A\nwMDBwcHCwsLDw8PExMTFxcXGxsbHx8fIyMjJycnKysrLy8vMzMzNzc3Ozs7Pz8/Q0NDR0dHS\n0tLT09PU1NTV1dXW1tbX19fY2NjZ2dna2trb29vc3Nzd3d3e3t7f39/g4ODh4eHi4uLj4+Pk\n5OTl5eXm5ubn5+fo6Ojp6enq6urr6+vs7Ozt7e3u7u7v7+/w8PDx8fHy8vLz8/P09PT19fX2\n9vb39/f4+Pj5+fn6+vr7+/v8/Pz9/f3+/v7///+V7mQAAAAACXBIWXMAABJ0AAASdAHeZh94\nAAAgAElEQVR4nO3df5xdZX3g8SNptEkErUWzKNoqaxakjSFITMOvWkKbSUChGJHwI2Kh2DUL\nFCssIAsWatOgLCTCKlVbLWxXCDYNxrbBgjW2soBI+LUshI3JzQwCYkogOf/ufb73uc+cZ77n\n3nme771xMjOfz+vl3HnOOfecb27O27nJJKEoiajnirEegGgiBCSiPgQkoj4EJKI+BCSiPgQk\noj4EJKI+BCSiPtQjpE8v7tKiRfLwu7Pnq13vmx0tB2YfNbw4fvZxXU46973dLjmyebMHuu0e\n6Lq32ZGzRzmg1aLZ81In6tJ7j6w58yLz6Y7MeqVGbd7s36vdbhxwweyFPQyjOnb28Z12DSwa\n7WdZd3zn00Wd2y9IH3+6y84XGrvdw6PFGWrX3CnRclexYHjxl8UXupx01hty5vtA8UK33Ttf\nGuX5c6YmXWZ7MZA6UZemH6427Wn81Hy6w17fyzCqxcW22u3P7TadbkVxXy/TjOz64quddu1s\n/Hv26b5YfDHlsF2Lw6dA6hqQQkDSAckHpPSApAOSD0jpAUkHJB+Q0gOSDkg+IKUHJB2QfEBK\nD0g6IPmAlB6QdEDyASk9IOmA5ANSekDSAckHpPSApAOSD0jpAUkHJB+Q0gOSDkg+IKUHJB2Q\nfEBKD0g6IPmAlB6QdEDyASk9IOmA5ANSekDSAckHpPSApAOSD0jpAUkHJB+Q0gOSDkg+IKUH\nJB2QfEBKD0g6IPmAlB6QdEDyASk9IOmA5ANSekDSAckHpPSApAOSD0jpAUkHJB+Q0gOSDkg+\nIKUHJB2QfEBKD0g6IPmAlB6QdEDyASk9IOmA5ANSekDSAckHpPSApNtHIS1tBaThgJQYkFoD\nAakEUi8BqTUQkEog9RKQWgMBqQRSLwGpNRCQSiD1EpBaAwGpBFIvAak1EJBKIPUSkFoDAakE\nUi8BqTUQkEog9RKQWgMBqQRSLwGpNRCQSiD1EpBaAwGpBFIvAak1EJBKIPUSkFoDAakEUi8B\nqTUQkEog9RKQWgMBqQRSLwGpNRCQSiD1EpBaAwGpBFIvAak1EJBKIPUSkFoDAakEUi8BqTUQ\nkEog9RKQWgMBqQRSLwGpNRCQSiD1EpBaAwGpBFIvAak1EJBKIPUSkFoDAakEUi8BqTUQkEog\n9RKQWgMBqQRSLwGpNVAV0rziCPfQ4aRASg1IrcYHpP+7p3MvNF51D5uLM9SuuVOi5cvFgubH\nGFKHk856Q5crqj5QPN9t989/Nsrz50xNusy2YiB1oi5NP1xt2t34qfl0h72+l2FUi4uf1G5/\n7lXT6VYU9/YyzchWFV/ptOvnjZ3Zp1tTrEk57OW+QVr+4FDndjTkYVNxmto1e0q03FbMa36M\nIXU46SEHdLmi6pji6W67B3eM8vzfnJp0mUeLhakTdWnaoXpbY7QJOzdrRg+z6E4sNtdub9hO\nd36xvodhVJ8tVnfaNdgYzD7dqmJVymHbFvULUn/e2i09rThwaYi3dsPx1i6xcf/WDkguIGUH\npCggSUDKDkhRQJKAlB2QooAkASk7IEUBSQJSdkCKApIEpOyAFAUkCUjZASkKSBKQsgNSFJAk\nIGUHpCggSUDKDkhRQJKAlB2QooAkASk7IEUBSQJSdkCK6hFS2w6QXECyN1kh/cprllYDkgtI\n9oAEpBCQ7AEJSCEg2QMSkEJAsjfpIHksQKoJSPaABKQQkOwBCUghINkDEpBCQLIHJCCFgGQP\nSEAKAckekIAUApI9IAEpBCR7QAJSCEj2gASkEJDsAQlIISDZAxKQQkCyByQghYBkD0hACgHJ\nHpCAFAKSPSABKQQke0ACUghI9oAEpBCQ7AEJSCEg2QMSkEJAsgckIIWAZA9IQAoByR6QgBQC\nkj0g7duQfmW/pV1mCQEpOyBFAUkCUnZAigKSBKTsgBQFJAlI2QEpCkgSkLIDUhSQJCBlB6So\ndEgRHCDVBSR7QAJSCEj2gASkEJDsAQlIISDZAxKQQkCyByQghYBkD0hACgHJHpCAFAKSPSAB\nKQQke0ACUghI9oAEpBCQ7AEJSCEg2QMSkEJAsgckIIWAZA9IQAoByR6QgBQCkj0gASkEJHtA\nAlIISPaABKQQkOwBCUghINkDEpBCQLIHJCCFgGQPSEAKAckekIAUApI9IAEpBCR7QAJSCEj2\ngASkEJDsAQlIISDZAxKQQkCyByQghYBkD0hACgHJHpCAFAKSPSABKQQke0ACUghI9oAEpBCQ\n7AEJSCEg2QPSuIDk63gYkLIDUhSQJCBlB6QoIElAyg5IUUCSgJQdkKKAJAEpOyBFAUkCUnZA\nigKSBKTsgBQFJAlI2QEpCkgSkLIDUhSQJCBlB6QoIElAyg5IUUCSgJQdkKKAJAEpOyBFAUkC\nUnZAigKSBKTsgBQFJAlI2QEpCkgSkLIDUhSQJCBlB6SoBEhLFxW/tnRkQKoJSPaABKQQkOwB\nCUghINkDEpBCQLIHJCCFgGQPSEAKAckekIAUApI9IAEpBCR7QAJSCEj2gASkEJDsAQlIISDZ\nAxKQQkCyByQghYBkD0hACgHJHpCAFAKSPSABKQQke0ACUghI9oAEpBCQ7AEJSCEg2QMSkEJA\nsgckIIWAZA9IQAoByR6QgBQCkj0gASkEJHtAAlIISPaABKQQkOwBCUghINkDEpBCQLK370Pa\nctES9/Di5878yJXbhh99QJKAlN1kg7Rx2UqBdNUlTz573QW7w6MPSBKQsptskDZsv89Baix+\novnV6KT724/t3UCSgJTdZINUlgLp3g/taX78xDfaj+2dQJKAlN0khbTuLPfppWvaj80P9xzb\n7LQHdnSp4T5YIXU45yEHdLviyI4pnuo6YGOU58eQOh62uTghZ6oOTTtUbxt1ws7NmtHDLLqF\nxSO1240Dnles62EY1TXFTR33GSZcWaxMOWzromxIZ3tIZwdImz7abOlDz3VusOE+WiF1OOkh\nB3S5ourY4pluuwcHR3l+DKnjYY8VC3Om6tC0Q/W2xg7z6WbN6GEW3YnFo7XbjQOeX9zdwzCq\nPy3WdNo11BjKPt2qYlXKYduzIX2v9ZbutvZjeydv7STe2mU3Sd/aDS5+rCyfX/JQ+7G9E0gS\nkLKbbJCGGuuXNBo7y2tWPLnligv3hEcfkCQgZTfZIC0fcN1RvrRy2elXD5Xh0QckCUjZTTZI\nowQkCUjZASkKSBKQsgNSFJAkIGUHpCggSUDKDkhRQJKAlB2QovYipA43L5BSA1IrIAEJSKkB\nCUghINkDEpBCQLIHJCCFgGQPSEAKAckekIAUApI9IAEpBCR7QAJSCEj2gASkEJDsAQlIISDZ\nAxKQQkCyByQghYBkD0hACgHJHpCAFAKSPSABKQQke0DKg7R06QFTu93RIwNS3wKSDkg+IKUH\nJB2QfEBKD0i68QQpvoeBlB6QWgEJSEBKDUhACgHJHpCAFAKSPSABKQQke0ACUghI9oAEpBCQ\n7AEJSCEg2QMSkEJAsgckIIWAZA9IQAoByR6QgBQCkj0gASkEJHtAAlIISPaABKQQkOwBCUgh\nINkDEpBCQLIHJCCFgGQPSEAKAckekIAUApI9IAEpBCR7QAJSCEj2gASkEJDsASkRUkggjXZj\nh8yQ2nMCqR2QdEDyASk9IOmA5ANSekDSAckHpPSApAOSD0jpAUkHJB+Q0gOSDkg+IKUHJB2Q\nfEBKD0i6iQJpdE9A6ltA0gHJB6T0gKQDkg9I6QFJByQfkNIDkg5IPiClByQdkHxASg9IOiD5\ngJQekHRA8gEpPSDpgOQDUnpA0gHJB6T0gKQDkg9I6QFJByQfkNIDkg5IPiClByQdkHxASg9I\nOiD5gJQekHRA8gEpPSDpgOQDUnpA0gHJB6T0gKQDkg9I6QFJByQfkNIDkg5IPiClByQdkHxA\nSg9IOiD5gJQekHRA8gEpPSDpJjSkaCOQ+haQdEDyASk9IOmA5ANSekDSAckHpPSApAOSD0jp\nAUk3MSHVbQRS3wKSDkg+IKUHJB2QfEBKD0g6IPmAlB6QdEDyASk9IOmA5ANSekDSAckHpPSA\npAOSD0jpAUk3wSDV1hoBSH0LSDog+YCUHpB0QPIBKT0g6YDkA1J6QNIByQek9ICkA5IPSOkB\nSQckH5DSA5IOSD4gpQckHZB8QEoPSDog+YCUHpB0QPIBKT0g6YDkA1J6QNIByQek9ICkA5IP\nSOkBSQckH5DSA5IOSD4gpQckHZB8QEoPSDog+YCUHpB0QPIBKT0g6YDkA1J6QNIByQek9ICk\nA5IPSOkBSQckH5DSA5JuMkBqNRNI/QpIOiD5gJQekHRA8gEpPSDpgOQDUnpA0gHJB6T0gKQD\nkg9I6QFJByQfkNIDkg5IPiClByQdkHxASg9IOiD5gJQekHRA8gEpPSDpgOQDUnpA0gHJB6T0\ngKQDkg9I6QFJByQfkNIDkq5/kJY/ONS5HQ33cYwhPd1lwKHBHZ321ELqeJ5Hi4XdLpPYtEP1\ntkbHCUdt1oweZtGdWGyu3d6wne78Yn0Pw6g+W6zutGuwMZh9ulXFqpTDti3qFyS+Ikl8RcqO\nr0hRQJKAlB2QooAkASk7IEUBSQJSdkCKApIEpOyAFAUkCUjZASkKSBKQsgNSFJAkIGUHpCgg\nSUDKDkhRQJKAlB2QooAkASk7IEUBSQJSdkCKApIEpOyAFAUkCUjZASkKSBKQsgNSFJAkIGUH\npCggSUDKDkhRQJKAlB2QooAkASk7IEUBSQJSdkCKApIEpOyAFAUkCUjZASkKSBKQsgNSFJAk\nIGUHpCggSUDKDkhRQJKAlB2QooAkASk7IEUBSQJSdkCKApIEpOyAFAUkCUjZTWhIcx5uPd4+\nK/WCQJKAlN2EhlRskodXrnxt6gWBJAEpuwkMqRjuvakXnMyQKhubkEazlhCQ7O1LkO5fVSxZ\n7vrY5c+kXhBIQAKSK3prd8KjuRcEEpCA5OJ37eoD0siApOsEaduZB+3X+kVS6gWBBCQguSIy\np/7S8WfKr5KWp14QSEACkiuC9KZv5l4QSEACkiuCNH177gWBBCQguSJIC/4h94JAAhKQXBGk\nHxx5b+YFgQQkILkiSPMPLqa/Q0q94L4P6ZRu93cmpBhLZQWk7CY0pAXHt0u9IJCABCTXhP+G\nLJD6FJB0QPIBKT0g6Tp+H6nd/qkXBBKQgOSKIC2Rjpx22AWpFwQSkIDkqntrt/XotakXBFIK\npG77KgHJ3j4Jqdw0J/WCQAISkFy1kLZOS70gkIAEJFcdpD3/7W2pFwQSkIDkiiD9hnTYrxYX\np14QSEACkqsG0uzjVr2cekEgAQlILr4hW197TiC1A5KuM6Qda9fcvK7rH5iOAxKQgOSKIO2+\naKr7BxtmXJt8QSABCUiuCNK1xcm3fGvt6hOKW1MvCCQgAckVQZp1Yevx4xPoX1oFUp8Ckq4T\npNdtaD3exTdk23MCqR2QdJ0gzbiz9fjN5J8CIAEJSK4I0m8dK99A2vk7x6ReEEhAApIrgnTX\na95+3lWfOfeg/e5OvSCQgAQkV/x9pP/1bvfb3++5K/mCQAISkFwj/2TDs9/f9JOMCwIJSEBy\nxZC2Xt/8sP3K+lesLiABCUiuCNIjM91/8/KpYuYTqReclJDi3NFAym5CQzrpXd93Dw+/64Op\nFwQSkIDkiiAd+D9aj6v5V4TacwKpHZB0nSBN88N8bXrqBYEEJCC5IkjvP+FV9/DC3PmpFwQS\nkIDkiiCte807L7ji8rMP3G9d6gWBBCQgueLf/l4/x31D9nC+IQskFZB0Xb4hu+OBH2X8BVkg\nAQlIEv9mQ33tOYHUDkg6IPmAlB6QdEDyASk9IOmA5ANSekDSAckHpPSApAOSD0jpAUkHJN/e\nhqSP7BaQ7AFJApILSPaAJAHJBSR7QJKA5AKSPSBJQHIByR6QJCC5gGQPSBKQXECyByRpPENy\nfbB4q/+s9iy1+3RAsgckCUguINkDkgQkF5DsAUkCkgtI9oAkAckFJHtAkoDkApI9IElAcgHJ\nHpAkILmAZA9I0sSB5IvPEm/sGJDsAUkCkgtI9oAkAckFJHtAkoDkApI9IElAcgHJHpAkILmA\nZA9IEpBcQLIHJAlILiDZA5IEJBeQ7AFJ2uuQOtzKQEoPSDog+YCUHpB0QPIBKT0g6YDkA1J6\nQNIByQek9ICkA5IPSOkBSQckH5DSA5IOSD4gpQckHZB8QEoPSDog+YCUHpB0QPIBKT0g6YDk\nA1J6QNIByQek9ICkA5IPSOkBSQckH5DSA5IOSD4gpQckHZB8QEoPSDog+YCUHpB0QPIBKT0g\n6YDkA1J6QNIByQek9ICkA5IPSOkBSQckH5DSA5IOSD4gpQckHZB8QEoPSDog+YCUHpB0Fkhb\nLlriHl783JkfuXL4FQWSBKTsJiukjctWCqSrLnny2esuCK8dkCQgZTdZIW3Yfp+D1Fj8RPOr\n0kn3tzcDSQJSdpMVUlkKpHs/tKf58RPfaG8EkgSk7CY5pHVnuU8vXdP8cM+xzU57YEeXGu7D\nvgRJDdjoNHp7TiOk+Cwdrj6yaYfWvIQdJxy1WTPMT61rYfFI7XbjgOcV63oYRnVNcVPHfYYJ\nVxYrUw7busgM6ewA6b7FzZY+ONS5HQ33cV+CNHLCwcFOs7fnNEKqrctLJU07VG9r7BjtWR2b\nNcP81LpOLDbXbjcOeH6xvodhVJ8tVnfaNdjo+LPcsVXFqpTDtpkhfa/11u629kbe2klpkEZ7\nhXlrZ2+8vbUbXPxYWT6/5KH2RiBJQMpuskIaaqxf0mjsLK9Z8eSWKy7c094MJAlI2U1WSMsH\nXHeUL61cdvrVQ2EzkCQgZTdZIXUISBKQsgNSFJAkIGUHpCggSUDKDkhRQJKAlB2QooAkASk7\nIEUBSQJSdkCKApIEpOyAFAUkCUjZASkKSBKQsgNSFJAkIGUHpCggSUDKDkhRQJKAlB2QooAk\nASk7IEUBSQJSdkCKApIEpOyAFAUkCUjZASkKSBKQsgNSFJAkIGUHpCggSUDKDkhRQJKAlB2Q\nooAkASk7IEUBSUqD1GGIEJDsAUkCkgtI9oAkAckFJHtAkoDkApI9IElAcgHJHpAkILmAZA9I\nEpBcQLIHJAlILiDZA5IEJBeQ7AFJApILSPaAJAHJBSR7QJKA5AKSPSBJQHIByR6QJCC5gGQP\nSBKQXECyByQJSC4g2QOSBCQXkOwBSQKSC0j2gCQByQUke0CSgOQCkj0gSUByAckekCQguYBk\nD0gSkFxAsgckCUguINkDkgQkF5DsAUkCkgtI9oAkAckFJHtAkoDkApI9IElAcgHJHpAkILmA\nZA9IEpBcQLIHJAlILiDZA5IEJBeQ7AFJApILSPaAJP2iII28o8cEUrj6yIBkD0gSkFxAsgck\nCUguINkDkgQkF5DsAUkCkgtI9oAkAckFJHtAkoDkApI9IElAcgHJHpAkILmAZA9IEpBcQLIH\nJAlILiDZA5IEJBeQ7AFJApILSPaAJAHJBSR7QJKA5AKSPSBJQHIByR6QJCC5gGQPSBKQXECy\nByQJSC4g2QOSBCQXkOwBSQKSC0j2gCQByQUke0CSgOQCkj0gSUByAckekCQguYBkD0gSkFxA\nsgckCUguINkDkgQkF5DsAUkCkgtI9oAkAckFJHtAkoDkApI9IElAcgHJHpAkILmAZA9IEpBc\nQLIHJAlILiDZA5IEJBeQ7AFJApILSPaAJAHJBSR7QJKA5AKSPSBJQHIByR6QJCC5gGQPSBKQ\nXECyByQJSC4g2QOSBCQXkOwBSQKSC0j2gCQByQUke0CSgOQCkj0gSUByAckekKTJCCkMEQKS\nPSBJQHIByR6QJCC5gGQPSBKQXECyByQJSC4g2QOSBCQXkOwBSQKSC0j2gCQByQUke0CSgOQC\nkj0gSUByAckekCQguYBkD0gSkFxAsgckCUguINkDkgQkF5DsAUkCkgtI9oAkAckFJHtAkoDk\nApI9IElAcgHJHpAkILmAZA9IEpBcQLIHJAlILiDZA5IEJBeQ7AFJApILSPaAJAHJBSR7QJKA\n5AKSPSBJQHIByR6QJCC5gGRvvEM654EdXWq4D/supPaEdbXn3OuQKtecdmjNS9hxwlGbNcP8\n1LoWFo/UbjcOeF6xrodhVNcUN3XcZ5hwZbEy5bCti/oFia9IEl+RsuMrUhSQJCBlB6QoIElA\nyg5IUUCSgJQdkKKAJAEpOyBFAUkCUnZAippAkMLG6gpIISDpgOTvYSClByQdkPw9DKT0gKQD\nkr+HgZQekHRA8vcwkNIDkg5I/h4GUnpA0gHJ38NASg9IOiD5exhI6QFJByR/DwMpPSDpgOTv\nYSClByQdkPw9DKT0gKQDkr+HgZQekHRA8vcwkNIDkg5I/h4GUnpA0gHJ38NASg9IOiD5exhI\n6QFJByR/DwMpPSDpgOTvYSClByQdkPw9DKT0gKQDkr+HgZQekHRA8vcwkNIDkg5I/h4GUnpA\n0gHJ38NASg9IOiD5exhI6QFJByR/DwMpPSDpgOTvYSClByQdkPw9DKT0gKQDkr+HgZQekHRA\n8vcwkNIDkg5I/h4GUnpA0gHJ38NASg9IOiD5exhI6QFJByR/DwMpPSDpgOTvYSClByQdkPw9\nDKT0gKQDkr+HgZQekHRA8vcwkNIDkg5I/h4GUnpA0gHJ38NASg9IOiD5exhI6QFJByR/DwMp\nPSDpgOTvYSClByQdkPw9DKT0gKQDkr+HgZQekHRA8vcwkNIDkg5I/h4GUnpA0gHJ38NASg9I\nOiBV89PWroAUApIOSNX8tLUrIIWApANSNT9t7QpIISDpgFTNT1u7AlIISDogVfPT1q6AFAKS\nDkjV/LS1KyCFgKQDUjU/be0KSCEg6YBUzU9buwJSCEg6IFXz09augBQCkg5I1fy0tSsghYCk\nA1I1P23tCkghIOmAVM1PW7sCUghIOiBV89PWroAUApIOSNX8tLUrIIWApANSNT9t7QpIISDp\ngFTNT1u7AlIISDogVfPT1q6AFAKSDkjV/LS1KyCFgKQDUjU/be0KSCEg6YBUzU9buwJSCEg6\nIFXz09augBQCkg5I1fy0tSsghYCkA1I1P23tCkghIOmAVM1PW7sCUghIOiBV89PWroAUApIO\nSNX8tLUrIIWApANSNT9t7QpIISDpgFTNT1u7AlIISDogVfPT1q6AFAKSDkjV/LS1KyCFgKQD\nUjU/be0KSCEg6YBUzU9buwJSCEg6IFXz09augBQCkg5I1fy0tSsghYCkA1I1P23tCkghIOmA\nVM1PW7sCUghIOiBV89PWroAUApIOSNX8tLUrIIWApANSNT9t7QpIISDpgFTNT1u7AlIISDog\nVfPT1q6AFAKSDkjV/LS1KyCFgKQDUjU/be0KSCEg6YBUzU9buwJSCEg6IFXz09augBQCkg5I\n1fy0tSsghYCkA1I1P23tCkghIOmAVM1PW7sCUghIOiBV89PWroAUApIOSNX8tLUrIIWApANS\nNT9t7QpIISDpgFTNT1u7AlIISDogVfPT1q6AFAKSDkjV/LS1KyCFgKQDUjU/be0KSCEg6YBU\nzU9buwJSCEg6IFXz09augBQCkg5I1fy0tSsghYCkA1I1P23tCkghIOmAVM1PW7sCUghIOiBV\n89PWroAUApIOSNX8tLUrIIWApANSNT9t7QpIISDpgFTNT1u7AlIISDogVfPT1q6AFAKSDkjV\n/LS1KyCFgKQDUjU/be0KSCEg6YBUzU9buwJSCEg6IFXz09augBQCkg5I1fy0tSsghYCkA1I1\nP23tCkghIOmAVM1PW7sCUghIOiBV89PWroAUApIOSNX8tLUrIIWApANSNT9t7QpIISDpgFTN\nT1u7AlIISDogVfPT1q6AFAKSrhdIWy5aUlkBSQJSdpMd0sZlK4GkAlJ2kx3Shu33AUkFpOwm\nO6SyBJIOSNkBqQ3pvsXNlj441LkdDfdxn4fkp61d7X1Ilddr2qH6NWzs6PICd2/WDPNT6zqx\n2Fy73Tjg+cX6HoZRfbZY3WnXYGMw+3SrilUph21b9AuENLSpOE3tmj0lHqiYN7y4sbi2y0kP\nOSDlR9jumOLpbrsHR3uJf3Nq0mUeLRamTtQlINkbx5CkhLd25aPFGWrX3CnRclexYHjxl8UX\nupx01hsSZgx9oHih2+6dL43y/DlTky6zvRhInahLvLWzN47f2klAkoCU3WSHNNRYv6TR2Nle\nAkkCUnaTHdLyAdcd7SWQJCBlN9khjQhIEpCyA1IUkCQgZQekKCBJQMoOSFFAkoCUHZCigCQB\nKTsgRQFJAlJ2QIoCkgSk7IAUBSQJSNkBKQpIEpCyA1IUkCQgZQekKCBJQMoOSFFAkoCUHZCi\ngCQBKTsgRQFJAlJ2QIoCkgSk7IAUBSQJSNkBKQpIEpCyA1IUkCQgZQekKCBJQMoOSFFAkoCU\nHZCigCQBKTsgRQFJAlJ2QIoCkgSk7IAUBSQJSNkBKQpIEpCyA1IUkCQgZQekKCBJQMoOSFFA\nkoCUHZCigCQBKTsgRX181Zc7t/qGL7mHP5v5u2rXEQdFy1tmzhte/OeZf9DlpIf/WpedqgUz\nV3fbffMXR3n+nLcmXeaGmcemTtSlg2erTV+64Sbz6X7j7b0Mozpu5hdqt9/4JdPpTpn5mV6m\nGdnHZ36y066bb7g5+3QXzLwg5bBb+gbpn/5nlz59zlfdw1dXXK12XbYiWt6+4vLhxY0rPt/l\npH9yUbdLjuyKFX/Tbfdtt43y/Mv+S9Jlvr7iM6kTdenCS9WmvznnYvPpLs16pUbtMyu+Xrv9\nG7ebTnfdii/3Ms3IPr/ixk67/vycv8g+3Q0rbkg67u/7Balrl8z5yV48+6Ro55zzxnqEcd9t\nc+78BVwFSPt0QOo9IBGQ+hCQCEh9aPxDIpo0AYmoDwGJqA8BiagP9RXSlouWqG0vfu7Mj1y5\nrXxgQFrbz8tNxLq8hGV518dO/sT3f/Ezjbe6vYZ/5O7CU/fCRfsJaeOylfpHcDmJkJAAAAXn\nSURBVNUlTz573QW7dzWa/ejUbn80j7q/hOW3l23a9s1zR/vjgZO+rq/h2Xc278PBvXDVfkLa\nsP0++REMXbvslE893trWWPxE8/8OTrpfFpf9VR+vNiHr+hKeu2EsRxs3dX0NT9m0l67a318j\ntX4EF137wstf+ejLsuXeD+1pfvzEN9znG5e/0terTcg6v4Q7Bjb80SkX/XgshxsndX4Ndw1c\n/8lzrt6yF665FyA9PjBUlns+vFG2rDvLfbx0TfPD7vPu7uvFJmadX8LNA3/yzAtrPmz/exWT\nps6v4U/P+PPNm68442f9v+ZegLSx9fsKt92zZMmSh9ed7bYLpI1nvdrXi03MOr+Emweab5Bf\nXfrtMR5wHNT1NizLn5+6vv/X3AuQvjfQ+nL60lNPPfXv32t9Tb2t+eHKNX291gSt80vYGHis\n+XjBbWM53fio623Y7A/3wq/V9wKkpwceaX7c2toyuLj5s//8kofK8mf+dxyoa51fwt3L7izL\nl39/41hONz7q/Bo+9fnmL9N3nvqd/l+zn5CGGuuXNBo7y0v/ePur3zrF/x7jNSue3HLFhc3/\nP7h/oP6vKlOlri/hbaf/sHH9sp1jO+G+X7fX8IWPrNy65eqz8//u+aj1E9JyeVN6Rzn0p6f9\n/h8/5De+tHLZ6Vc3f9lX/sNifs9u1Lq+hLtvPePkT/GduNHq+ho+cdlpH71qb/ylBP6IEFEf\nAhJRHwISUR8CElEfAhJRHwISUR8CElEfAhJRHwLSmPS+/ziRLkNAGqNWXl1d/TD7Z2H0Z8gR\n8WWyz0Hp8WLuA12f/bMw+jP6cQSlx4s5Jrn3XAt+69+O2//AD28rTyiKYk5Z/uNv7z9t9i3N\nnfMX3Pm2eWW5/ujXv+VU9zcnwo73HrVh7rQ3nv3T9jOGC4f8v4+9/XVv+eCP/RGty2yc+8sH\nXbvrkoNef/wTzSP+eu60/ef8damuSj0FpDHJ3eHHHzz37m23TzmzfHRJsenh8ttTjr5z/XnF\nn5XlcYe/+4a15frX/M5Xb3nnf9ha2THvwCP+ufGVqSf7Zww3fMhRM2/+ztfe8+aXWkfIZd52\n7L8+c3Lx21du+acDfq8sv16cvHbtwmKtuir1FJDGJLnDi+82Pzv+oLJc7n4WZr/L/ftAi/ff\n2dzxt83Pjvj1V8ryX167qrJjfuH+MtLy4unWM4YLhzxffKr5yeNXP9s6onWZ+8vynuL9zeXp\nM8ry6uNeLsvnf+l0dVXqKSCNSXKHT3efnblf65beVnxyZ7Obiu+Xx792V1nuKP6wdWhlx/wZ\n7u953lp8awSk4UN2vekd3279pyiHITX1lI8XFzc/Xhz+Q6BvW6CuSj0FpDFJ7vB3uM/c7ez+\n98PC97fyRap8sLiidWhlx/x3ug1riy+NgFQ55Lu/XrzpQ197pQrJXeb/FNc0P15SPFc+/18P\nO2DKlGK+uir1FJDGpDpI59wnNVo7flRc3jq0sqMF6ZvFrQpSOKR8dcPF/6k44ucdIR095dMb\nH3jwoPnqqtRTQBqTNKTB4sz2TtnxQiH/8M1T2ys75k9z/w7T6mLdCEiVQ6T/Xny5E6THinOb\nn7zyy/PVVamngDQmxZA+VjTfix35hueay1svfcXveM+BzV/Q/Lj5Bm94x/zmr47K8qTXDbWe\nMVw45AenuX8Y4/HiutYRGtLDxZWl+xbSUeqq1FNAGpNiSJcXV95e/uPUw2/9+8umntXesXa/\n9/3VmkPevLWyY/7Bh9x498XFMv+M4cIhW/c//Ja7v/7+Ax5vHaEh7Tr4rXd896Jjjtn/Oz8b\ncVXqKSCNSTGkZ2ZPba7v+cD+Uw+59pX2jvKuo6a/+eRHy8qO+e/+wdHT33jui+1nDBcO+d8n\nv3nqQSf/mz+i5tdIm+ZNf8sfPH/nr75x84irUk8BaRw1nz+Dus8GpHEUkPbdgDSOAtK+G5DG\nURGkvytCN47ZRNQOSOO1Fx8MDY31LAQkon4EJKI+BCSiPgQkoj4EJKI+BCSiPgQkoj70/wFB\n5hB3yY+uMAAAAABJRU5ErkJggg==",
520 "output_type": "display_data"
524 "ggplot(df, aes(x=intercept_estimate)) + geom_histogram(bins=100) +\n",
525 " geom_vline(xintercept=dh$xbr) + \n",
526 " scale_x_log10() + scale_y_log10() +\n",
528 "ggsave(paste0('plot_dhist_', unique(df$op), '.png'))"
539 "codemirror_mode": "r",
540 "file_extension": ".r",
541 "mimetype": "text/x-r-source",
543 "pygments_lexer": "r",