BAYES: Hồi quy nhị thức âm
Sử dụng ngôn ngữ STAN
Nhóm tác giả BAV
07 Tháng 11 2017
Tác giả:
Đinh Tiến Tài (1,2,3) Lê Đông Nhật Nam (1,3,4) Trần Trung Dũng (3) Lê Ngọc Khả Nhi (1,4,5)
Ghi chú: Đóng góp của các thành viên:
- Bài giảng lý thuyết
- Mô hình Poisson và NBI cổ điển
- Soạn và hoàn thiện STAN codes cho mô hình NBI Bayes
- Khai thác mô hình và Suy diễn Bayes
- Đồ họa thống kê
1 Giới thiệu
Thân chào các bạn, đây là bài thực hành thứ 7 trong dự án Bayes for Vietnam với sự góp mặt của toàn thể coreteam của project. Mục tiêu của nhóm BAV là phổ cập về phương pháp thống kê theo trường phái Bayes nhằm thay thế hoàn toàn những công cụ truyền thống. Đối tượng của chúng tôi là các bạn bác sĩ và sinh viên y khoa.
Đây cũng là một bài giảng với tham vọng cao nhất,vì chúng tôi sẽ đề cập đồng thời nhiều vấn đề tương đối phức tạp, bao gồm: Phân tích biến số đếm (count data) bằng phân phối Nhị thức âm, mixed model có chứa random effect, và Hồi quy tuyến tính đa biến (những bài giảng trước kia chỉ dừng lại ở phân tích đơn biến).
Chúng tôi giả định rằng các bạn đã bắt đầu quen với cấu trúc ngôn ngữ STAN, quy trình chuyển giả thuyết nghiên cứu thành mô hình Bayes, khai thác phân phố hậu định .Những phần này sẽ được giản lược
2 Trường hợp minh họa
library(tidyverse)
library(ggridges)
my_theme <- function(base_size = NULL, base_family = "sans") {
theme_bw(base_size = base_size, base_family = base_family) +
theme(
panel.grid.major = element_line(color = "gray"),
panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "#fff68f" ),
strip.background = element_rect(fill = "#510166", color = NA, size =0.5),
strip.text = element_text(face = "bold", size = 8, color = "white"),
legend.position = "bottom",
legend.justification = "center",
legend.background = element_blank(),
legend.margin = margin(0.5,0.5,0.5,0.5)
)
}Bộ dữ liệu này được thu thập từ một nghiên cứu trên 59 bệnh nhân mắc bệnh động kinh. Các bệnh nhân được theo dõi số lần bị động kinh trong vòng 8 tuần (base). Sau đó bệnh nhân được phân chia ngẫu nhiên vào nhóm điều trị (treatment- dùng thuốc Progabide ) và nhóm giả dược (placebo); Số lần bị động kinh được ghi nhận mỗi 2 tuần trong 8 tuần tiếp theo. Như vậy mỗi bệnh nhân sẽ được khảo sát ở 4 thời điểm (period) và dạng dữ liệu này được gọi là dữ liệu tái đo lường (repeated measurement data). Dataset đã được sử dụng bởi nhiều tác giả như Thall ,Vail [1990], Breslow và Clayton[1993], Lee và Nelder [1996, 2000] để minh họa cho mô hình hồi quy hỗn hợp và count data.
Đây là một thử nghiệm lâm sàng tiêu biểu, trong đó câu hỏi nghiên cứu là khảo sát tác động của một loại thuốc X làm thay đổi (được kì vọng là tích cực) một kết cục xấu của một bệnh lý. Thiết kế nghiên cứu thông dụng là phân chia ngẫu nhiên và đồng nhất một mẫu n bệnh nhân vào 2 phân nhóm (nhánh) trị liệu, một cho dùng thuốc X, một sử dụng giả dược (placebo). Kết cục (outcome) Y sẽ được theo dõi kéo dài qua nhiều thời điểm trong suốt quá trình điều trị. Một mô hình hồi quy đa biến trong đó phân nhóm điều trị (X vs placebo) là biến số chính, kèm theo những hiệp biến như thời gian, tuổi, tình trạng cơ bản (baseline)… và yếu tố ngẫu nhiên cá thể/thứ bậc khác. Suy diễn thống kê được thực hiện trên tham số hồi quy cho X.
Trong trường hợp này, kết cục lâm sàng Y là số cơn động kinh được ghi nhận tính đến mỗi mốc thời gian (k=4), yếu tố can thiệp có 2 levels là thuốc progabide so với placebo. Giả thuyết nghiên cứu là thuốc Progabide cải thiện có ý nghĩa (làm giảm số cơn động kinh) so với Placebo. Bài toán của chúng ta còn xét thêm 2 hiệp biến số tiềm ẩn khác đó là : Tình trạng cơn động kinh (lbase : logarit của tần suất động kinh trước khi phân nhóm ) và Thời gian (visit : lưu ý, trong thí nghiệm gốc, thời gian được đo như một biến liên tục, tính bằng tháng, và có 4 thời điểm khảo sát, tuy nhiên chúng tôi hoán chuyển biến số này thành những trọng số tương phản : contrast weights) để đơn giản hóa mô hình. Việc xét thêm 2 biến số này nhằm để hiệu chỉnh sai lệch (nếu có) do độ nặng cá thể lên suy diễn thống kê về hiệu quả điều trị, cũng như hiệu ứng phụ của thời gian.
dat=read.csv("https://raw.github.com/vincentarelbundock/Rdatasets/master/csv/MASS/epil.csv")
names(dat)[2] <- "sei"
dat$visit <- (2*dat$period-5)/10
dat$period=as.factor(dat$period)
dat$subject=as.factor(dat$subject)
dat%>%head()%>% knitr:: kable()| X | sei | trt | base | age | V4 | subject | period | lbase | lage | visit |
|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 5 | placebo | 11 | 31 | 0 | 1 | 1 | -0.7563538 | 0.1142037 | -0.3 |
| 2 | 3 | placebo | 11 | 31 | 0 | 1 | 2 | -0.7563538 | 0.1142037 | -0.1 |
| 3 | 3 | placebo | 11 | 31 | 0 | 1 | 3 | -0.7563538 | 0.1142037 | 0.1 |
| 4 | 3 | placebo | 11 | 31 | 1 | 1 | 4 | -0.7563538 | 0.1142037 | 0.3 |
| 5 | 3 | placebo | 11 | 30 | 0 | 2 | 1 | -0.7563538 | 0.0814139 | -0.3 |
| 6 | 5 | placebo | 11 | 30 | 0 | 2 | 2 | -0.7563538 | 0.0814139 | -0.1 |
Nhưng trước khi đi vào bài toán, chúng tôi muốn bàn một chút về lý thuyết xác suất và mô hình cho biến số đếm (count data)
Thống kê mô tả cho thấy Outcome (Seizure count) có phân phối không bình thường : Cụ thể, nó có trung bình nhỏ hơn Sd, skewness rất cao.
psych::describe(dat$sei)## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 236 8.25 12.35 4 5.83 4.45 0 102 102 4.16 22.06 0.8sigma=sd(dat$sei)
mu=mean(dat$sei)
realp=dat%>%ggplot(aes(x=sei))+
geom_density(alpha=.8,fill="#c507ff")+
geom_histogram(aes(y=..density..,fill=..density..),colour="black",alpha=0.7,show.legend = F,binwidth = 1)+
theme_bw()+scale_fill_gradient(low="#ff9c07",high="#ff072f")+
ggtitle("Real data")
realp
Nguyên nhân vì đây là một số đếm (count data).
Trong y học lâm sàng, nhiều đại lượng được khảo sát như một số đếm / biến số rời rạc (count data, discrete variables). Chúng thường mang ý nghĩa tần suất (số lần) một biến cố/hiện tượng được ghi nhận trong một khoảng không gian và thời gian xác định, Thí dụ :
Hiện tượng sinh lý : nhịp tim, nhịp hô hấp, nhu động ruột, số cơn gò tử cung,
Biến cố bệnh lý : Số cơn hen kịch phát, số lần ngưng thở khi ngủ, số cơn động kinh, số trường hợp tử vong, số lần nhập viện cấp cứu…
Kết quả của phép đếm : Số lượng tế bào bạch cầu, số kí sinh trùng, số con, số đơn vị cơ quan bị tổn thương…
Hành vi : số lần thăm khám bác sĩ, mức tiêu thụ thuốc lá/bia rượu, số lần dùng thuốc cắt cơn hen
Hoặc một thang đo rời rạc : điểm Glasgow, mức độ đau,
Thời gian cũng có thể được xét như biến số đếm, khi chúng là số nguyên, thí dụ: Số ngày nằm viện, độ dài giấc ngủ ,…
Đặc điểm cần lưu ý với count data là phân bố dữ liệu thường hay bị lệch phải, phương sai tăng khi giá trị trung bình tăng.
3 Biến số đếm: từ đại lượng đến biến ngẫu nhiên
Tiếp theo, chúng tôi muốn phân biệt giữa đại lượng và biến số trong nghiên cứu, vì có sự khác nhau giữa 2 khái niệm này
Xác định được bản chất, ý nghĩa và cách đo đếm những đại lượng là một bước quan trọng khi đặt giả thuyết và thiết kế nghiên cứu, đây là trách nhiệm của người bác sĩ. Nhưng bước tiếp theo, quan trọng không kém, đó là chuyển từ đại lượng sang biến số –và đây là công việc của chuyên viên thống kê.
Một cách lý tưởng, người bác sĩ có đầy đủ cả kiến thức thống kê và y học, sẽ cùng lúc đảm nhận cả 2 vai trò để tự mình đi từ quan sát đến giả thuyết, từ giả thuyết đến thí nghiệm, đo/đếm đại lượng và từ đại lượng này đến biến số thống kê.
Các mô hình thống kê là những con robot vô tri, chúng hoàn toàn không “nhìn” thấy đại lượng, không biết đến ý nghĩa lâm sàng. Chúng chỉ biết tới biến số ngẫu nhiên. Chúng không đo lường, mà chỉ dùng hàm toán học để ước lượng giá trị của biến số cần xét. Mỗi giá trị như vậy lại có xác suất hiện diện khác nhau. Nhà thống kê không có kiến thức y học chuyên môn sẽ không quan tâm đến ý nghĩa sinh lý, bệnh học nữa, nhưng chỉ cố gắng mô tả đặc tính phân phối của biến ngẫu nhiên bằng quy luật phù hợp (hàm mật độ xác suất).Các hàm này được định nghĩa bằng các tham số và giả định. Từ phân phối, họ có thể ước lượng tham số/giá trị bằng mô hình, và suy diễn thống kê.
Tới đây, các bạn sinh viên và bác sĩ cũng có thể nhận ra vấn đề : Những biến số đếm, rời rạc không thể được mô tả bằng cùng một phân phối Gaussian (chuẩn, bình thường ) như những biến số liên tục. Như vậy các phương pháp xử lý số liệu có giả định phân phối chuẩn (t-test, ANOVA, OLS regression, etc) khi áp dụng với dữ liệu count data sẽ không còn phù hợp nữa (trừ khi dữ liệu đã được chuyển đổi). Khi bạn đã giả định sai về quy luật phân phối, mô hình của bạn đã có bias và nhiều nguy cơ là suy diễn thống kê dựa trên mô hình này cũng sai lầm nốt.
Điều đáng tiếc đó là chương trình Thống kê ở đại học Y chỉ quan tâm đến những biến liên tục, nhưng chưa trang bị cho các sinh viên Y khoa đủ kiến thức để xử trí những biến rời rạc và biến số đếm.
Một mâu thuẫn khác, đó là khi thi nội trú, cao học, người ta có dạy lý thuyết xác suất về những quy luật phân phối như Binomial (Nhị thức ), Bernoulli, Poisson, nhưng cách dạy thuần túy lý thuyết và toán học không cho phép vận dụng những quy luật phân phối này vào thực tiễn lâm sàng và nghiên cứu – sinh viên không có khả năng liên hệ sự vật, hiện tượng với quy luật phân phối, và không có khả năng chuyển từ đại lượng thành biến số ngẫu nhiên.
Một sai lầm chết người khác, đó là việc dạy cho sinh viên tiếp cận vấn đề bằng tư duy so sánh, và các công cụ lạc hậu như bảng chéo, kiểm định phi tham số… nhưng không dạy về mô hình hồi quy, nhất là mô hình tuyến tính tổng quát (GLM).
Với 2 lỗ hổng kiến thức cả về Xác suất ứng dụng, Mô hình GLM, sinh viên gần như hoàn toàn bất lực trước những giả thuyết nghiên cứu với biến số đếm. Sinh viên chỉ có 2 lựa chọn, làm liều và sai (dùng kiểm định t so sánh, vốn dựa trên giả định phân phối chuẩn), hoặc tìm nơi trú ẩn (các test phi tham số).
Trong dự án Bayes analysis for Vietnam, nhóm chúng tôi nhắm đến cả 2 mục tiêu : Mục tiêu chính là khuyến khích các bạn đi theo trường phái Bayes, nhưng một mục tiêu khác quan trọng không kém, đó là bù đắp cho các bạn đồng nghiệp những thiếu hụt về kiến thức và kỹ năng thống kê ứng dụng. Mỗi bài thực hành là một case study hoàn chỉnh dẫn dắt các bạn từ câu hỏi, giả thuyết đến biến số, mô hình và suy diễn kết quả.
4 Từ phân phối nhị thức đến Nhị thức âm
Phân phối nhị thức và Bernoulli là gốc rễ của mọi quy luật phân phối dành cho biến rời rạc.
Phân phối nhị thức mô tả biến số ngẫu nhiên x như Số lần thành công trên n lần thử nghiệm lặp lại độc lập (n xác định), biết rằng xác suất thành công là p:
\[P(x|n,p)=\binom{n}{x}p^{x}(1-p)^{n-x}\]
Phân phối nhị thức âm :
Dạng chính tắc của phân phối nhị thức âm mô tả số lần thử nghiệm Bernoulli độc lập, lặp lại (biến ngẫu nhiên x) để có thể đạt được một số lần thành công xác định (lần thành công thứ r), biết rằng xác suất thành công là p
\[P(x|p,r)= \binom{x-1}{r-1}p^{r}(1-p)^{x-r}\] Với x và r đều là số nguyên, điều kiện x >= r, và p là một xác suất trong khoảng (0,1).
Dạng thứ hai của phân phối nhị thức âm là một biến thể, trong đó biến ngẫu nhiên x là số lần thử nghiệm thất bại trước khi lần thành công thứ r xảy ra. \[P(x|p,r)= \binom{x+r-1)}{x}p^{r}(1-p)^{x}\] Với x là số nguyên dương hoặc=0 và không có ràng buộc về điều kiện giữa r và x. Như vậy, chính dạng biến thể này mới thực sự hữu dụng để mô tả biến số đếm vì x có thể nhận bất cứ giá trị số nguyên nào không âm, và r cũng không bị ràng buộc điều kiện nào cả.
Tên gọi “nhị thức âm có vẻ bí ẩn ? Thực ra nó phát xuất từ trường hợp ta áp dụng phân phối nhị thức với số mũ âm :
\[\binom{x+r-1}{x}\ p^{r}(1-p)^{x} = \binom{-r}{x}\ p^{r}(p-1)^{x}\]
Thay vì biểu diễn phân phối nhị thức âm bằng 2 tham số r và p, ta có thể sử dụng khái niệm khác là : tham số vị trí trung tâm (Mu) và phương sai (sigma^2), quan hệ hoán chuyển giữa Mu, Sigma^2 và r,p như sau :
\[\mu = \frac {r(1-p)}{p}\]
\[Var = \sigma^2 = \frac{r(1-p)}{p^2} = \mu + \frac{1}{r}\mu^2\] Dạng này có tên gọi là Negative binomial type 2 hoặc quadratic negative binomial. Số 2 chính là bậc mũ 2 cho Mu trong công thức tính variance.
5 Biện luận về phân phối cho bài toán minh họa
Trở lại bài toán của chúng ta:
Những lập luận sau đây ủng hộ cho việc chọn họ phân phối phù hợp cho mô hình :
Về mặt lâm sàng : Outcome của chúng ta có bản chất là tần suất phát sinh một biến cố lâm sàng (cơn động kinh), được khảo sát trong một khoảng thời gian xác định, trong một không gian xác định (cơ thể mỗi bệnh nhân).
Một số ý khác có thể được phát biểu, như: Mỗi cơn động kinh là độc lập (thực vậy !), ngẫu nhiên (giả định) và có xác suất hằng định theo thời gian (giả đinh). Những đặc tính này thỏa đa số điều kiện và tính chất của một phân phối Poisson dành cho biến số đếm, nhưng không phù hợp với phân phối Gaussian. Do đó, phân phối Poisson là giả định ban đầu.
Khi các bạn muốn xem xét ảnh hưởng của các biến độc lập x1,x2,..xn lên biến phụ thuộc dạng số đếm y thì chúng ta thường hay nghĩ ngay đến mô hình hồi quy Poisson. Mô hình Poisson có một giả định (assumption) rất quan trọng là: giá trị mean và variance bằng nhau. Tuy nhiên dữ liệu thực tế mà chúng ta thu thập ít khi thõa mãn giả định này. Khi variance > mean chúng ta gọi đó là over-dispersion (tạm dịch: quá phân tán), ngược lại, variance < mean thì gọi là under-dispersion.
Một biến ngẫu nhiên y với giá trị kì vọng μ tuân theo phân bố poisson sẽ có hàm probability mass function (hàm khối xác suất) như sau:
\[P\left ( Y = y \right ) = \frac{e^{-\mu }\mu^y}{y!}\] Có thể nhận ra rằng phân phối poisson chỉ có 1 tham số là mean (mu). Trong hồi quy poisson, mu được mô tả bằng một hàm số mũ của các biến độc lập tuyến tính (linear predictor). Link function được sử dụng trong hồi quy poisson là log link.
\[\mu = e^\eta\] \[\eta = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + ...+ \beta_n x_n\]
Đầu tiên chúng ta sẽ dùng mô hình poisson không có random effect để kiểm tra độ phân tán của dữ liệu (dispersion).
glm.pois <- glm(sei ~ lbase + trt + visit , family = "poisson", data = dat) # fitting mô hình
summary(glm.pois)##
## Call:
## glm(formula = sei ~ lbase + trt + visit, family = "poisson",
## data = dat)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -5.4185 -1.4334 -0.2444 0.7166 10.9926
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.76538 0.03823 46.174 < 2e-16 ***
## lbase 1.17646 0.03086 38.121 < 2e-16 ***
## trtprogabide -0.10325 0.04532 -2.278 0.02272 *
## visit -0.29598 0.10148 -2.917 0.00354 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2517.83 on 235 degrees of freedom
## Residual deviance: 973.81 on 232 degrees of freedom
## AIC: 1747.7
##
## Number of Fisher Scoring iterations: 5library(COUNT)
P__disp(glm.pois)## pearson.chi2 dispersion
## 1120.845606 4.831231Hàm P__disp cung cấp cho chúng ta Pearson Chi2 statistic và Pearson dispersion statistic. Có thể tính thủ công Pearson dispersion statistic bằng cách lấy tổng bình phương của phần dư Pearson chia cho bậc tự do của mô hình (N quan sát - số biến trong mô hình bao gồm intercept).
poi.df <- df.residual(glm.pois) # 232
Pear.disp.sta <- sum(resid(glm.pois, type = "pearson")^2)/poi.df
Pear.disp.sta## [1] 4.831231Giá trị Pearson dispersion statistic > 1, có nghĩa là dữ liệu gặp vấn đề overdispersion. Xin có một lưu ý với các bạn là khi dữ liệu bị overdispersion, một giải pháp thay thế cho mô hình poisson là mô hình quasi-poisson.
glm.quasi <- glm(sei ~ lbase + trt + visit , family = "quasipoisson", data = dat)
summary(glm.quasi) %>% coef()## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.7653831 0.08403771 21.007035 1.362195e-55
## lbase 1.1764626 0.06783359 17.343363 9.064533e-44
## trtprogabide -0.1032456 0.09961567 -1.036439 3.010764e-01
## visit -0.2959814 0.22304668 -1.326993 1.858156e-01Hệ số hồi quy của mô hình glm.poi và glm.quasi là giống nhau, standard error (SE) của các hệ số hồi quy ở mô hình glm.pois lớn hơn so với SE ở mô hình glm.pois.
Chúng ta có thể tính toán sự sai khác SE giữa 2 mô hình dựa vào Pearson dispersion statistic. Khi nhân SE của các hệ số hồi quy ở mô hình poisson với căn bậc hai của Pearson dispersion statistic chúng ta sẽ có được SE tương ứng trong mô hình quasi-poisson.
Những giải pháp thay thế khác có thể là : phân phối Poisson Inversed Gaussian (PIG), hoặc Negative binomial type II (Nhị thức âm, NBII).
Chúng ta có thể kiểm chứng về tính phù hợp của phân phối negative binomial II bằng cách : Đối chiếu histogram cho dữ liệu thực tế và histogram cho một dữ liệu mô phỏng theo quy luật NBI2, với cùng giá trị tham số mu và sigma. Kết quả cho thấy 2 đồ thị gần như tương đương :
library(gamlss.dist)
rNBII(236L, mu = mu, sigma = sigma)%>%as_tibble()%>%
ggplot()+
geom_density_ridges(aes(x=value,y=1,group=1),stat = "binline",binwidth=1,
draw_baseline = FALSE,alpha=.8,fill="gold")+
geom_density_ridges(aes(x=dat$sei,y=1,group=1),stat = "binline",binwidth=1,
draw_baseline = FALSE,alpha=.5,fill="#c507ff",col="#c507ff")+
theme_bw()+scale_y_continuous("Density")+
ggtitle("Simulated (yellow) vs Real data (violet)")
Nhận xét : Việc lựa chọn họ phân phối có thể dựa vào một hoặc cả 2 lập luận : 1) Bản chất/Ý nghĩa lâm sàng và 2) Toán thống kê.
Với (1), bạn đang nhìn đại lượng, với đầy đủ ý nghĩa, đơn vị, thang đo… theo quan điểm y học, với (2), outcome không còn là đại lượng nữa mà chỉ là một biến số ngẫu nhiên, giá trị của nó là 1 con số vô nghĩa chịu chi phối bởi các hàm mật độ xác suất.
Trên thực tế hiếm khi ta có thể liên hệ một cách tường minh giữa ý nghĩa lâm sàng và định nghĩa xác suất, thí dụ tỉ lệ điều trị thành công trên tổng số bệnh nhân vừa mang ý nghĩa lâm sàng, vừa đúng theo định nghĩa của một phân phối Binomial ; nhưng trong trường hợp này, rất khó để hình dung số cơn động kinh là số lần thử nghiệm Bernouilli thất bại cho đến khi đạt lần thành công thứ 1/alpha - theo định nghĩa phân phối nhị thức âm type 2 ; tương tự, phân phối Poisson có thể được lựa chọn chỉ vì ta không có thông tin nào khác, chứ trên thực tế không phải bất cứ đại lượng nào cũng thỏa tất cả giả định về tính ngẫu nhiên, độc lập và bất biến của phân phối Poisson. Thí dụ : mùa xuân là thời điểm bệnh nhân hen dị ứng với phấn hoa khởi phát nhiều cơn hen hơn so với mùa hè, xác suất của cơn hen đã bị thay đổi chứ không còn hằng định nữa.
Như vậy, phân phối nhị thức âm được chọn cho nghiên cứu này như một phương tiện, để điều chỉnh tính quá phân tán của số liệu, chứ không phải như một quy luật.
5.1 Mô hình NBII : Từ REML đến Bayes
Cho bài toán hồi quy với mục tiêu ước tính giá trị biến kết quả y (thay cho x) có phân phối NBI với tham số trung bình Mu và kiểu hình 1/alpha (alpha thay cho 1/r, theo định nghĩa về số lần thành công thứ r). Phân phối Negative binomial còn có thể được trình bày dưới hình thức
\[P(y|\mu,\alpha) = \binom{y+1/\alpha -1}{1/\alpha-1}( \frac{1}{1+\alpha \mu})^{1/\alpha} (\frac{\alpha \mu}{1+\alpha \mu})^{y}\]
Hay hình thức hỗn hợp hàm Poisson và Gamma :
\[P(y|\mu,\alpha) = \frac {\Gamma(y+1/\alpha)}{\Gamma(y+1)\Gamma(1/\alpha)}( \frac{1}{1+\alpha \mu})^{1/\alpha} (\frac{\alpha \mu}{1+\alpha \mu})^{y}\] Điều kiện duy nhất ở đây là Mu và Alpha > 0. Như đã nói ở trên, 1/alpha = r không bắt buộc phải là số nguyên .
Mô hình negative binomial type 2 còn được gọi mô hình poisson-gamma distribution mixture model. Với biến y tuân theo phân phối NB2 thì chúng ta cần lưu ý những điểm sau:
Giả định y tuân theo phân phối poisson với mean = variance = Mu
Mu sẽ tuân theo phân phối gamma với Mean = Mu và
\[variance = \mu + \frac{\mu^2}{v}\]
Ở đây v được gọi là shape parameter (tham số kiểu hình) trong phân phối gamma, nó tương ứng với tham số r trong phân phối NBI chính tắc (còn gọi là indirect dispersion parameter trong phân phối NB); ta thay r hay v bằng alpha là direct dispersion parameter.
Như vậy đối với mô hình NB, có 2 tham số cần được ước tính là Mu (trung bình) và dispersion parameter (v hoặc alpha). Lưu ý rằng khi sử dụng R, thì v được gọi bằng tên thay thế là Theta ở trong phần kết quả của mô hình. Một số phần mềm khác như Stata, SAS không tính toán theta hay v mà sẽ tính alpha và gọi đó là dispersion parameter trong phần kết quả.
Khi NB2 model được mô hình hóa dưới dạng kết hợp của gamma-poisson thì variance của mô hình NB sẽ bằng tổng variance của phân phối poisson và phân phối gamma.
Tham số Mu có thể được ước lượng bằng mô hình GLM thông qua 1 link function là hàm logarit tự nhiên (LN) :
\[ln(\mu) \sim \beta_{0} + \sum_{i=1}^{K}\beta_{i} x_{i}\]
Với beta0 là intercept và beta(i) là tham số hồi quy cho biến số độc lập x(i) trong model matrix X gồm k biến số
Theo trường phái ước tính hợp lý cực đại (Maximum likelihood), ta có thể ước lượng tham số alpha và beta (cho mu) thông qua hàm likelihood :
\[L(\beta|\alpha ,\beta )= \prod_{i=1}^{n}P(y_{i}) = \prod_{i=1}^{n} \frac {\Gamma(y_{i}+1/\alpha)}{\Gamma(y_{i}+1)\Gamma(1/\alpha)}( \frac{1}{1+\alpha e^{x_{i}\beta}})^{1/\alpha} (\frac{\alpha e^{x_{i}\beta}}{1+\alpha e^{x_{i}\beta}})^{y_{i}}\]
Giá trị alpha và beta hợp lý nhất tương ứng với L đạt cực đại.
Mô hình mixed NB với random effect theo trường phái frequentist có thể được fit với package gamlss hay lme4. Để đơn giản về mặt cú pháp, chúng ta sẽ dùng package lme4. Với bộ dữ liệu epil, chúng ta sẽ xem xét ảnh hưởng của trt, lbase và visit lên số lần bị động kinh với random intercept cho từng bệnh nhân.
require(lme4)
nb.admb <- glmer.nb(sei ~ trt+lbase+visit + (1|subject), data = dat )
summary(nb.admb)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: Negative Binomial(7.4538) ( log )
## Formula: sei ~ trt + lbase + visit + (1 | subject)
## Data: dat
##
## AIC BIC logLik deviance df.resid
## 1264.6 1285.4 -626.3 1252.6 230
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.9341 -0.5748 -0.0558 0.4678 3.5659
##
## Random effects:
## Groups Name Variance Std.Dev.
## subject (Intercept) 0.2391 0.4889
## Number of obs: 236, groups: subject, 59
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.8055 0.1078 16.742 <2e-16 ***
## trtprogabide -0.3332 0.1510 -2.207 0.0273 *
## lbase 1.0092 0.1006 10.032 <2e-16 ***
## visit -0.2714 0.1656 -1.639 0.1013
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) trtprg lbase
## trtprogabid -0.708
## lbase -0.028 -0.064
## visit 0.010 0.002 -0.005Khi dựng mô hình hồi quy nhị thức âm theo Bayes, mục tiêu của chúng ta là ước tính phân phối hậu định cho tham số beta (cho phép ước tính Mu) và tham số alpha là kiểu hình. Mô hình này cũng có thể chứa random effect. Prior và likelihood của mô hình Bayes NBI được minh họa trong hình sau:
Mô hình mixed NB với random intercept chúng ta cần ước tính có dạng:
\[\mu_{ij} = e^\left ( \beta_{0} + \beta_{trt} x_{trt} + \beta_{lbase} x_{lbase} + \beta_{visit} x_{visit} + sub_i + \varepsilon_{ij} \right )\]
6 Lựa chọn mô hình NBI2 trong STAN
Bây giờ chúng ta sẽ sử dụng Stan để để fit mô hình thông qua package rstan.
- Bước 1: Viết Stan code.
Ở những bài thực hành trước các bạn đã làm quen với cách viết stan code theo từng block, nên chúng tôi xin phép không nhắc lại lần nữa. Ở đây chúng tôi sẽ ghi chú trực tiếp trong từng block để giải thích ý nghĩa của từng dòng lệnh.
stan_NBI.stri ="
data {
int<lower=0> N; //no of obs
int<lower=0> K; //no of parameters (intercept+predictors)
int<lower=0> Nsub; //no of subs
int<lower=0> sub[N]; // index for subject
int<lower=0> y[N]; //dependent variable
matrix[N,K] X; //model matrix
vector[K] bpri;
matrix[K,K] Bpri;
vector[Nsub] apri;
matrix[Nsub, Nsub] Apri;
}
parameters {
vector[K] beta; //beta[1] is population intercept
vector[Nsub] re_itcept; // random intercept
real<lower=0> sigma_re; //SD random intercept
real<lower=0> alpha; //direct dispersion parameter
}
transformed parameters {
vector[N] log_lik;
vector[N] eta; //linear predictor
vector[N] mu; //estimated mean of dependent variable
eta = X*beta;
for(i in 1:N) {
mu[i] = exp(eta[i] + re_itcept[sub[i]]);
log_lik[i] = neg_binomial_2_lpmf(y[i]|mu[i],1/alpha);
}
}
model {
//prior
alpha ~ cauchy(0,25);
sigma_re ~ cauchy(0,5);
beta ~ multi_normal(bpri, Bpri); //mean = 0; sd = 1
re_itcept ~ multi_normal(apri, Apri*sigma_re); //mean = 0
//likelihood
y ~ neg_binomial_2(mu, 1/alpha); //theta = 1/alpha
}
"- Bước 2: Lựa chọn mô hình bằng vòng lặp WAIC
Trước khi thực sự dựng mô hình, chúng tôi sẽ giới thiệu về phương pháp lựa chọn mô hình (biến số) dựa vào tiêu chí WAIC (Watanabe-Ikaike information criteria) và LOOIC (Leave-one-out information criteria). Cơ chế đằng sau những trị số này sẽ được giải thích rõ hơn trong những bài tiếp theo.
Tại thời điểm này, chúng ta tạm chấp nhận quy tắc : WAIC và LOOIC được tính từ Log-likelihood của mô hình Bayes, và mô hình phù hợp nhất với dữ liệu sẽ có WAIC và LOOIC thấp nhất. bằng cách so sánh WAIC và LOOIC giữa các phiên bản model khác nhau, ta có thể chọn được model hợp lý nhất. Cách làm này giống như quy trình Step-wise dựa vào AIC hay BIC trong trường phái Maximum Likelihood mà bạn đã biết.
Chúng tôi viết 1 vòng lặp, lần lượt dựng 7 mô hình NBI Bayes cho tất cả tổ hợp có thể giữa 3 biến số: Treatment, Logbase và Visit, và tính WAIC, LOOIC cho mỗi mô hình.
library(loo)## This is loo version 1.1.0library(rstan)
library(parallel)
com <- c(combn(2:4, 1, simplify=FALSE),
combn(2:4, 2, simplify=FALSE),
combn(2:4, 3, simplify=FALSE))
trt=dat$trt
lbase=dat$lbase
visit=dat$visit
Xmat = model.matrix(~ trt+lbase+visit)
output=data.frame(
WAIC=rep(NA,7),
SEWAIC=rep(NA,7),
LOOIC=rep(NA,7),
SELOOIC=rep(NA,7),
ModelID=c(1:7)
)
Nsub = nlevels(as.factor(dat$subject)) # no of subjects
N=nrow(dat)
Pbar <- winProgressBar(title = "WAIC loop", min = 0,
max = length(com), width = 350)
for(m in 1:length(com)) {
X= Xmat[,c(1, unlist(com[[m]]))]
K= ncol(X)
tempdat <- list(y = dat$sei, # response variable
N=N, # no of obs
X=X, # model matrix of predictors and intercept
K=K, # no of predictors + intercept
Nsub = Nsub, # no of subs ( 59 = no of random effect VALUEs)
sub = dat$subject%>%as.numeric(), # Subject ID
bpri = rep(0, K), # To define prior for beta
Bpri = diag(1, K), # Prior matrix of beta
apri = rep(0, Nsub), # To define prior for intercept
Apri = diag(1, Nsub )) # Prior matrix of intercept
set.seed(123)
tempmod <- stan(model_code = stan_NBI.stri, data = tempdat,
iter = 1500, warmup = 500)
w1=tempmod%>%extract_log_lik(.,parameter_name ="log_lik")%>%waic()
l1=tempmod%>%extract_log_lik(.,parameter_name ="log_lik")%>%loo()
output$WAIC[m]<-w1$waic
output$SEWAIC[m]<-w1$se_waic
output$LOOIC[m]<-l1$looic
output$SELOOIC[m]<-l1$se_looic
setWinProgressBar(Pbar, m, title=paste( round(m/length(com)*100, 0),
"% Hoàn thành"))
}## In file included from C:/Users/Admin/Documents/R/win-library/3.4/BH/include/boost/config.hpp:39:0,
## from C:/Users/Admin/Documents/R/win-library/3.4/BH/include/boost/math/tools/config.hpp:13,
## from C:/Users/Admin/Documents/R/win-library/3.4/StanHeaders/include/stan/math/rev/core/var.hpp:7,
## from C:/Users/Admin/Documents/R/win-library/3.4/StanHeaders/include/stan/math/rev/core/gevv_vvv_vari.hpp:5,
## from C:/Users/Admin/Documents/R/win-library/3.4/StanHeaders/include/stan/math/rev/core.hpp:12,
## from C:/Users/Admin/Documents/R/win-library/3.4/StanHeaders/include/stan/math/rev/mat.hpp:4,
## from C:/Users/Admin/Documents/R/win-library/3.4/StanHeaders/include/stan/math.hpp:4,
## from C:/Users/Admin/Documents/R/win-library/3.4/StanHeaders/include/src/stan/model/model_header.hpp:4,
## from file20a8330a4d06.cpp:8:
## C:/Users/Admin/Documents/R/win-library/3.4/BH/include/boost/config/compiler/gcc.hpp:186:0: warning: "BOOST_NO_CXX11_RVALUE_REFERENCES" redefined
## # define BOOST_NO_CXX11_RVALUE_REFERENCES
## ^
## <command-line>:0:0: note: this is the location of the previous definition
## cc1plus.exe: warning: unrecognized command line option "-Wno-ignored-attributes"
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 11.619 seconds (Warm-up)
## 13.366 seconds (Sampling)
## 24.985 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 11.715 seconds (Warm-up)
## 21.832 seconds (Sampling)
## 33.547 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 11.062 seconds (Warm-up)
## 26.736 seconds (Sampling)
## 37.798 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 11.219 seconds (Warm-up)
## 26.077 seconds (Sampling)
## 37.296 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 7.743 seconds (Warm-up)
## 11.114 seconds (Sampling)
## 18.857 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0.01 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 100 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 8.067 seconds (Warm-up)
## 11.136 seconds (Sampling)
## 19.203 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 7.811 seconds (Warm-up)
## 11.119 seconds (Sampling)
## 18.93 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 7.677 seconds (Warm-up)
## 11.058 seconds (Sampling)
## 18.735 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 9.265 seconds (Warm-up)
## 11.207 seconds (Sampling)
## 20.472 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 8.493 seconds (Warm-up)
## 11.195 seconds (Sampling)
## 19.688 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 8.015 seconds (Warm-up)
## 11.208 seconds (Sampling)
## 19.223 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 8.136 seconds (Warm-up)
## 11.208 seconds (Sampling)
## 19.344 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 9.621 seconds (Warm-up)
## 11.379 seconds (Sampling)
## 21 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 9.796 seconds (Warm-up)
## 11.125 seconds (Sampling)
## 20.921 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 10.184 seconds (Warm-up)
## 12.231 seconds (Sampling)
## 22.415 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 10.005 seconds (Warm-up)
## 12.633 seconds (Sampling)
## 22.638 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 14.654 seconds (Warm-up)
## 24.133 seconds (Sampling)
## 38.787 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 13.355 seconds (Warm-up)
## 24.452 seconds (Sampling)
## 37.807 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 13.315 seconds (Warm-up)
## 27.658 seconds (Sampling)
## 40.973 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 12.269 seconds (Warm-up)
## 25.732 seconds (Sampling)
## 38.001 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 7.633 seconds (Warm-up)
## 12.656 seconds (Sampling)
## 20.289 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 8.523 seconds (Warm-up)
## 11.777 seconds (Sampling)
## 20.3 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 8.411 seconds (Warm-up)
## 12.658 seconds (Sampling)
## 21.069 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 9.853 seconds (Warm-up)
## 13.673 seconds (Sampling)
## 23.526 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0.002 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 20 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 12.302 seconds (Warm-up)
## 13.781 seconds (Sampling)
## 26.083 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 9.966 seconds (Warm-up)
## 11.844 seconds (Sampling)
## 21.81 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 11.128 seconds (Warm-up)
## 11.554 seconds (Sampling)
## 22.682 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0.001 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 1500 [ 0%] (Warmup)
## Iteration: 150 / 1500 [ 10%] (Warmup)
## Iteration: 300 / 1500 [ 20%] (Warmup)
## Iteration: 450 / 1500 [ 30%] (Warmup)
## Iteration: 501 / 1500 [ 33%] (Sampling)
## Iteration: 650 / 1500 [ 43%] (Sampling)
## Iteration: 800 / 1500 [ 53%] (Sampling)
## Iteration: 950 / 1500 [ 63%] (Sampling)
## Iteration: 1100 / 1500 [ 73%] (Sampling)
## Iteration: 1250 / 1500 [ 83%] (Sampling)
## Iteration: 1400 / 1500 [ 93%] (Sampling)
## Iteration: 1500 / 1500 [100%] (Sampling)
##
## Elapsed Time: 10.224 seconds (Warm-up)
## 11.775 seconds (Sampling)
## 21.999 seconds (Total)library(viridis)
output%>%ggplot(aes(x=reorder(ModelID,WAIC),y=WAIC))+
geom_point(shape=21,aes(fill=factor(ModelID),
size=WAIC,col=factor(ModelID),
stroke=SEWAIC/1.5),
show.legend = F,alpha=0.5)+
geom_point(shape=21,aes(fill=factor(ModelID),
size=WAIC),
show.legend = F)+
theme_bw()+scale_x_discrete("Model version")+
scale_colour_viridis(option="C",discrete=T)+
scale_fill_viridis(option="C",discrete = T)
output%>%ggplot(aes(x=reorder(ModelID,LOOIC),y=LOOIC))+
geom_point(shape=21,aes(fill=factor(ModelID),
size=LOOIC,col=factor(ModelID),
stroke=SELOOIC/1.5),
show.legend = F,alpha=0.5)+
geom_point(shape=21,aes(fill=factor(ModelID),
size=LOOIC),
show.legend = F)+
theme_bw()+scale_x_discrete("Model version")+
scale_colour_viridis(option="C",discrete=T)+
scale_fill_viridis(option="C",discrete = T)
Kết quả cho thấy: Mô hình phù hợp nhất là mô hình số 7, nội dung của nó là:
colnames(X)## [1] "(Intercept)" "trtprogabide" "lbase" "visit"Như vậy ta có thể chắc chắn trong mô hình này có 3 biến số như giả thuyết ban đầu
Trước khi dựng mô hình, chúng ta có thể thăm dò dữ liệu 1 chút:
- Tương phản giữa 2 phân nhóm trị liệu tại 4 thời điểm (visit khác nhau:)
den.rid <- dat %>% ggplot(., aes(x = sei, y = trt, fill= trt)) +
geom_density_ridges(stat = "binline", scale = 1, binwidth=2, draw_baseline = FALSE) +
labs(x="Count", y = "") + scale_fill_discrete(name = "Treatment") +
coord_flip() +
my_theme(10) +
geom_vline(xintercept = median(dat$base),linetype=2,col="blue")+
theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1)) +
theme(axis.text = element_text(size = 12, color = "black"),
axis.title = element_text(size = 15, color = "black")) +
facet_wrap(~ period, ncol=4)
den.rid 
Ta có thể thấy rằng tình trạng bệnh lý trước khi điều trị (log baseline)có ảnh hưởng ít nhiều đến hiệu quả điều trị của progabide:
g.smooth <- dat%>%ggplot(aes(x=trt,y=sqrt(sei), group=lbase))+
geom_point(alpha=0.8,aes(color=lbase))+
geom_smooth(size=1,aes(color=lbase),se=F,method="lm",show.legend = F)+
my_theme(12)+scale_color_viridis(option="A", direction = 1)+
theme(axis.text = element_text(size = 11, color = "black"),
axis.title = element_text(size = 15, color = "black")) +
facet_wrap(~ period,ncol=4)
g.smooth 
Ta có thể khảo sát hiệu ứng của Thời gian lên Logarit của biến số kết quả, có vẻ như hiệu ứng này rất yếu:
g.trt <- dat %>% ggplot(.,aes(x=period, y=log(sei + 0.01), group = trt)) +
geom_point(aes(color = trt), size = 2 ) +
geom_point(shape =1, size = 2, color = "black") +
geom_smooth(aes(color= trt, fill = trt), size = 1, alpha= 0.3, se = T, method="lm") +
my_theme(12) +
theme(axis.text = element_text(size = 12, color = "black"),
axis.title = element_text(size = 15, color = "black"))
g.trt 
Cuối cùng, ta có thể khảo sát khuynh hướng đáp ứng điều trị ở mỗi cá thể bệnh nhân, bây giờ thì bạn có thể hiểu tại saota dùng mô hình với random effect
dat%>%ggplot(aes(x=period,y=sei,group=1))+
geom_path(aes(color=trt),alpha=0.9,size=1)+
geom_point(aes(color=trt,shape=trt),alpha=0.9,size=2)+
my_theme(5)+
facet_wrap(~subject,scales="free",ncol=10)
7 Mô hình NBII Bayes chính thức
Xmat <- model.matrix(~dat$trt+dat$lbase+dat$visit)
Nsub <- nlevels(dat$subject) # no of subjects
K = ncol(Xmat)
data.NBI = list(y = dat$sei, # response variable
N = nrow(dat), # no of obs
X = Xmat, # model matrix of predictors and intercept
K = K, # no of predictors + intercept
Nsub = Nsub, # no of subs ( 59 = no of random effect VALUEs)
sub = dat$subject %>% as.numeric(.), # Subject ID
bpri = rep(0, K), # To define prior for beta
Bpri = diag(1, K), # Prior matrix of beta
apri = rep(0, Nsub), # To define prior for intercept
Apri = diag(1, Nsub )) # Prior matrix of intercept fit.stan <- stan(model_code = stan_NBI.stri,
data = data.NBI,
chains=4,
iter=3000,
warmup=1000,
thin=2)##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 1).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 3000 [ 0%] (Warmup)
## Iteration: 300 / 3000 [ 10%] (Warmup)
## Iteration: 600 / 3000 [ 20%] (Warmup)
## Iteration: 900 / 3000 [ 30%] (Warmup)
## Iteration: 1001 / 3000 [ 33%] (Sampling)
## Iteration: 1300 / 3000 [ 43%] (Sampling)
## Iteration: 1600 / 3000 [ 53%] (Sampling)
## Iteration: 1900 / 3000 [ 63%] (Sampling)
## Iteration: 2200 / 3000 [ 73%] (Sampling)
## Iteration: 2500 / 3000 [ 83%] (Sampling)
## Iteration: 2800 / 3000 [ 93%] (Sampling)
## Iteration: 3000 / 3000 [100%] (Sampling)
##
## Elapsed Time: 18.506 seconds (Warm-up)
## 23.931 seconds (Sampling)
## 42.437 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 2).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 3000 [ 0%] (Warmup)
## Iteration: 300 / 3000 [ 10%] (Warmup)
## Iteration: 600 / 3000 [ 20%] (Warmup)
## Iteration: 900 / 3000 [ 30%] (Warmup)
## Iteration: 1001 / 3000 [ 33%] (Sampling)
## Iteration: 1300 / 3000 [ 43%] (Sampling)
## Iteration: 1600 / 3000 [ 53%] (Sampling)
## Iteration: 1900 / 3000 [ 63%] (Sampling)
## Iteration: 2200 / 3000 [ 73%] (Sampling)
## Iteration: 2500 / 3000 [ 83%] (Sampling)
## Iteration: 2800 / 3000 [ 93%] (Sampling)
## Iteration: 3000 / 3000 [100%] (Sampling)
##
## Elapsed Time: 17.089 seconds (Warm-up)
## 22.948 seconds (Sampling)
## 40.037 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 3).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 3000 [ 0%] (Warmup)
## Iteration: 300 / 3000 [ 10%] (Warmup)
## Iteration: 600 / 3000 [ 20%] (Warmup)
## Iteration: 900 / 3000 [ 30%] (Warmup)
## Iteration: 1001 / 3000 [ 33%] (Sampling)
## Iteration: 1300 / 3000 [ 43%] (Sampling)
## Iteration: 1600 / 3000 [ 53%] (Sampling)
## Iteration: 1900 / 3000 [ 63%] (Sampling)
## Iteration: 2200 / 3000 [ 73%] (Sampling)
## Iteration: 2500 / 3000 [ 83%] (Sampling)
## Iteration: 2800 / 3000 [ 93%] (Sampling)
## Iteration: 3000 / 3000 [100%] (Sampling)
##
## Elapsed Time: 17.702 seconds (Warm-up)
## 22.783 seconds (Sampling)
## 40.485 seconds (Total)
##
##
## SAMPLING FOR MODEL 'f3d81cce00e79713520bd97b4c1e119a' NOW (CHAIN 4).
##
## Gradient evaluation took 0 seconds
## 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Adjust your expectations accordingly!
##
##
## Iteration: 1 / 3000 [ 0%] (Warmup)
## Iteration: 300 / 3000 [ 10%] (Warmup)
## Iteration: 600 / 3000 [ 20%] (Warmup)
## Iteration: 900 / 3000 [ 30%] (Warmup)
## Iteration: 1001 / 3000 [ 33%] (Sampling)
## Iteration: 1300 / 3000 [ 43%] (Sampling)
## Iteration: 1600 / 3000 [ 53%] (Sampling)
## Iteration: 1900 / 3000 [ 63%] (Sampling)
## Iteration: 2200 / 3000 [ 73%] (Sampling)
## Iteration: 2500 / 3000 [ 83%] (Sampling)
## Iteration: 2800 / 3000 [ 93%] (Sampling)
## Iteration: 3000 / 3000 [100%] (Sampling)
##
## Elapsed Time: 16.989 seconds (Warm-up)
## 23.811 seconds (Sampling)
## 40.8 seconds (Total)Kết quả của mô hình Bayes NBI2 được tóm tắt như sau:
print(fit.stan, digits=3, pars=c("beta","sigma_re","alpha"))## Inference for Stan model: f3d81cce00e79713520bd97b4c1e119a.
## 4 chains, each with iter=3000; warmup=1000; thin=2;
## post-warmup draws per chain=1000, total post-warmup draws=4000.
##
## mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff
## beta[1] 1.791 0.003 0.116 1.559 1.713 1.792 1.869 2.015 1614
## beta[2] -0.309 0.004 0.160 -0.628 -0.415 -0.311 -0.200 0.002 1700
## beta[3] 1.001 0.002 0.108 0.786 0.927 1.000 1.073 1.205 1885
## beta[4] -0.268 0.003 0.168 -0.597 -0.377 -0.267 -0.156 0.057 4000
## sigma_re 0.286 0.002 0.078 0.166 0.230 0.276 0.331 0.467 2697
## alpha 0.145 0.001 0.033 0.089 0.121 0.141 0.166 0.218 2932
## Rhat
## beta[1] 1.002
## beta[2] 1.001
## beta[3] 1.001
## beta[4] 1.000
## sigma_re 1.000
## alpha 1.000
##
## Samples were drawn using NUTS(diag_e) at Tue Nov 07 23:00:58 2017.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).Ta chỉ quan tâm đến phân phối hậu định của 3 hiệu ứng chính,dưới dạng Incidence rate ratio (IRR = exp(beta))
postdf=as.data.frame(fit.stan)%>%as_tibble()
betadf=postdf[,c(1:4)]
names(betadf)=c("Intercept","Progabide","LogBase","Visit")
betadf=betadf%>%mutate(.,Iteration=as.numeric(rep(c(1:1000),4)),
Chain=as.factor(rep(c(1:4),each=1000)))
head(betadf)%>%.[,c(1:4)]%>%knitr::kable()| Intercept | Progabide | LogBase | Visit |
|---|---|---|---|
| 1.772972 | -0.2198768 | 0.9795190 | -0.2581374 |
| 1.941806 | -0.3589941 | 1.0442256 | -0.3001051 |
| 1.862598 | -0.5882802 | 0.9337627 | -0.2383386 |
| 1.914448 | -0.7154537 | 1.0115780 | -0.1133664 |
| 1.942097 | -0.4211530 | 0.9456440 | -0.3032589 |
| 1.906347 | -0.4570913 | 0.8811814 | 0.0534731 |
p1=betadf%>%gather(Intercept:Visit,key="Step",value="fix")%>%
ggplot(aes(x=exp(fix),fill=Chain,col=Chain))+
geom_density(alpha=0.3,show.legend = F)+
facet_wrap(~Step,ncol=1,scales = "free_y")+
scale_x_continuous("IRR")+
theme_bw(6)
p2=betadf%>%gather(Intercept:Visit,key="Step",value="fix")%>%
ggplot(aes(y=exp(fix),x=Iteration,col=Chain))+
geom_path(alpha=0.5,show.legend = F)+
facet_wrap(~Step,ncol=1,scales = "free")+
scale_y_continuous("IRR")+
theme_bw(6)
gridExtra::grid.arrange(p2,p1,ncol=2)
Nếu chỉ xét Incidence rate change cho yếu tố điều trị, ta có thể diễn giải:Thuốc progabid làm giảm tần suất cơn động kinh khoảng 26.72 % (4.36 % đến 43.7% ) so với nhóm Placebo.
Hmisc::describe(exp(betadf$Progabide)-1)## exp(betadf$Progabide) - 1
## n missing distinct Info Mean Gmd .05 .10
## 4000 0 4000 1 -0.2564 0.134 -0.43730 -0.40309
## .25 .50 .75 .90 .95
## -0.33991 -0.26720 -0.18105 -0.09964 -0.04358
##
## lowest : -0.5879168 -0.5842844 -0.5775129 -0.5679095 -0.5621748
## highest: 0.2723041 0.2743226 0.2746589 0.2774731 0.2935617Tương quan giữa Incidence rate change của Visit và Progabide được biểu diễn như sau:
betadf$pseudoGroup=factor(rep(c(1:20),e=nrow(betadf)/20))
betadf%>%gather(Visit,Progabide,key="Factors",value="fix")%>%
ggplot(aes(y=Factors,
x=100*(exp(fix)-1)
))+
geom_density_ridges_gradient(aes(fill = ..x..),
show.legend = F,
scale=0.8,
gradient_lwd = 0.5)+
geom_density_ridges(aes(col=pseudoGroup),
scale=0.9,alpha=0.01,
show.legend = F)+
geom_vline(xintercept=c(0,-10,10),col=c("red3","blue3","black"),linetype=2)+
scale_x_continuous("Incidence rate change (%)",
expand = c(0.01, 0),
breaks=c(-70,-60,-50,-40,-30,-20,-10,0,10,20,30,40,50))+
scale_y_discrete(expand = c(0.01, 0))+
scale_fill_viridis(option = "A",begin=0.2,end=1)+
scale_colour_viridis(option = "A",discrete = T)+
theme_bw()
Phân tích ROPE theo J.Jruschke với ngưỡng trên và dưới lần lượt là -5% và +5% , cùng phân tích CompVal với ngưỡng là 10%
HDIF= function( sampleVec,credMass=0.975 ) {
sortedPts = sort( sampleVec )
ciIdxInc = ceiling( credMass * length( sortedPts ) )
nCIs = length( sortedPts ) - ciIdxInc
ciWidth = rep( 0 , nCIs )
for ( i in 1:nCIs ) {
ciWidth[ i ] = sortedPts[ i + ciIdxInc ] - sortedPts[ i ]
}
HDImin = sortedPts[ which.min( ciWidth ) ]
HDImax = sortedPts[ which.min( ciWidth ) + ciIdxInc ]
HDIlim = c( HDImin , HDImax )
return( HDIlim )
}
SUMK=function(paramSampleVec,compVal=NULL , ROPE=NULL , credMass=0.975) {
meanParam = mean( paramSampleVec )
medianParam = median( paramSampleVec )
dres = density( paramSampleVec )
modeParam = dres$x[which.max(dres$y)]
hdiLim = HDIF( paramSampleVec , credMass=credMass )
if ( !is.null(compVal) ) {
pcgtCompVal = ( 100 * sum( paramSampleVec > compVal )
/ length( paramSampleVec ) )
} else {
compVal=NA
pcgtCompVal=NA
}
if ( !is.null(ROPE) ) {
pcltRope = ( 100 * sum( paramSampleVec < ROPE[1] )
/ length( paramSampleVec ) )
pcgtRope = ( 100 * sum( paramSampleVec > ROPE[2] )
/ length( paramSampleVec ) )
pcinRope = 100-(pcltRope+pcgtRope)
} else {
ROPE = c(NA,NA)
pcltRope=NA
pcgtRope=NA
pcinRope=NA
}
return( c( Mean=meanParam , Median=medianParam , Mode=modeParam ,
HDIlevel=credMass , LL=hdiLim[1] , UL=hdiLim[2] ,
CompVal=compVal , PcntGtCompVal=pcgtCompVal ,
ROPElow=ROPE[1] , ROPEhigh=ROPE[2] ,
PcntLtROPE=pcltRope , PcntInROPE=pcinRope , PcntGtROPE=pcgtRope ) )
}
summaryKruschke=function(MCMC,compVal=NULL, rope=NULL,credMass=NULL){
summaryInfo = NULL
summaryInfo = cbind(summaryInfo, "Estimated"= SUMK(MCMC,
compVal=compVal,
ROPE=rope,credMass=credMass))
return(summaryInfo)
}
IRC=100*(1-exp(betadf$Progabide))%>%as_data_frame()
summaryKruschke(MCMC = IRC$value,
compVal=10,
rope=c(-5,5),
credMass=0.975)%>%as.data.frame()## Estimated
## Mean 25.637268
## Median 26.719580
## Mode 29.149447
## HDIlevel 0.975000
## LL -2.373055
## UL 50.106432
## CompVal 10.000000
## PcntGtCompVal 89.950000
## ROPElow -5.000000
## ROPEhigh 5.000000
## PcntLtROPE 1.325000
## PcntInROPE 4.075000
## PcntGtROPE 94.600000betadf$IRC=100*(1-exp(betadf$Progabide))Kết quả cho thấy: 94.6% mật độ phân phối hậu định nằm ngoài ngưỡng thay đổi 5%, chỉ có 4.075% mật độ rơi vào trong khoảng +/-5%; 89.95% mật độ nằm ngoài ngưỡng 10%
Để chắc chắn, ta làm thêm 1 phân tích Bayes Factor cho 21 ngưỡng thay đổi của IR từ 0 đén-20%
library(brms)
betadf$IRC=100*(1-exp(betadf$Progabide))
threshold=rep(NA,21)
BayesFactor=rep(NA,21)
thres=c(0:20)
for(i in (1:21)){
thr=thres[i]
threshold[i]=thr
hyp=paste("IRC>",thr,sep="")
bf=brms::hypothesis(betadf,hyp,alpha=0.05)
BayesFactor[i]=bf$hypothesis$Evid.Ratio
}
bfdf=cbind(threshold,BayesFactor)%>%as_tibble()
bfdf%>%ggplot(aes(x=threshold,y=BayesFactor,fill=BayesFactor))+
geom_path()+
geom_point(show.legend = F,size=5,shape=21,col="black")+
geom_text(aes(label=round(BayesFactor,2)),col="black",show.legend = F,angle = 60,nudge_y=sqrt(bfdf$BayesFactor)+1,nudge_x=0,size=4)+
theme_bw()+scale_x_continuous(breaks=c(0:20))+
geom_hline(yintercept = c(10,30),linetype=2,col="blue")+
scale_fill_viridis(option="D",direction=-1)
Nếu lấy ngưỡng khả tín của BF=30, ta có thể kết luận là tuy Progabide có cải thiện triệu chứng động kinh, nhưng hiệu ứng này khá nhỏ, Với H1 là IR thay đổi > 10% thì BF chỉ có 8.95; Ta chỉ có thể tin cậy hiệu ứng điều trị từ 1-5% mà thôi.
8 Diễn đạt văn bản khoa học
Phương pháp thống kê :
Do kết quả (số cơn động kinh) là một biến số đếm và bị phân tán, nó được khảo sát theo quy luật phân phối nhị thức âm type II (NBII). Hiệu ứng của thuốc Progabide so với Placebo lên kết cục lâm sàng được khảo sát bằng một mô hình hồi quy hỗn hợp (Mixed model) theo phương pháp Bayes. Mô hình này có hiệu chỉnh cho tình trạng bệnh tại khởi điểm nghiên cứu (lbase), yếu tố thời gian (4 trọng số tương phản) và hiệu ứng ngẫu nhiên cá thể. Tính hợp lý của mô hình được đánh giá bằng tiêu chí WAIC. Suy diễn thống kê dựa vào mật độ phân phối hậu định và tỉ trọng chứng cứ (Bayes factor) với giả thuyết H1 là thuốc Progabide giảm tần suất động kinh ít nhất là 5% so với Placebo.
Kết quả:
Phân phối hậu định của sự thay đổi tần suất phát sinh cơn động kinh (IR change,%) cho thấy thuốc Progabid cải thiện trung bình 26.72% số cơn động kinh so với nhóm Placebo. Đến 94.6% mật độ phân phối hậu định nằm ngoài ngưỡng thay đổi 5%, chỉ có 4.075% mật độ rơi vào trong khoảng +/-5%; 89.95% mật độ nằm ngoài ngưỡng 10. Tuy nhiên Bayes Factor cho thấy khả năng xác tín về hiệu quả này là khá thấp (BF=17.52 cho ngưỡng 5% và 8.95 cho ngưỡng 10%).
9 Tổng kết
Bài thực hành số 7 đến đây là chấm dứt. Trong bài này chúng tôi đã truyền tải một số thông điệp chính như sau:
Khác với thí nghiệm y học lâm sàng, đối tượng của phân tích thống kê không phải là những đại lượng sinh lý bệnh, mà là những biến số ngẫu nhiên đại diện cho các đại lượng này. Việc lựa chọn quy luật phân phối phù hợp để chuyển từ đại lượng sang biến số yêu cầu người bác sĩ phải vận dụng lý thuyết xác suất.
Những biến số đếm không thể được khảo sát bằng phương pháp truyền thống dựa trên giả định phân phối Gaussian. Một số giải pháp đã được giới thiệu, bao gồm mô hình Poisson, Quasi-Poisson và Nhị thức âm type 2.
Phân tích Bayes trong STAN cho phép dựng mô hình GLM với phân phối uyển chuyển, như trong trường hợp này là phân phối NBII. Suy diễn Bayes cũng mang lại nhiều thông tin hơn các kiểm định cổ điển.
Nhóm BAV xin chân thành cảm ơn sự ủng hộ của các bạn. Chúng tôi sẽ gặp lại các bạn trong một bài khác.
LS0tDQp0aXRsZTogIkJBWUVTOiBI4buTaSBxdXkgbmjhu4sgdGjhu6ljIMOibSINCnN1YnRpdGxlOiAiU+G7rSBk4bulbmcgbmfDtG4gbmfhu68gU1RBTiINCmF1dGhvcjogIk5ow7NtIHTDoWMgZ2nhuqMgQkFWIg0KZGF0ZTogIjA3IFRow6FuZyAxMSAyMDE3Ig0Kb3V0cHV0Og0KICBodG1sX2RvY3VtZW50OiANCiAgICBjb2RlX2Rvd25sb2FkOiB0cnVlDQogICAgY29kZV9mb2xkaW5nOiBoaWRlDQogICAgbnVtYmVyX3NlY3Rpb25zOiB5ZXMNCiAgICB0aGVtZTogImRlZmF1bHQiDQogICAgdG9jOiBUUlVFDQogICAgdG9jX2Zsb2F0OiBUUlVFDQotLS0NCg0KYGBge3Igc2V0dXAsaW5jbHVkZT1GQUxTRX0NCmtuaXRyOjpvcHRzX2NodW5rJHNldChlY2hvID0gVFJVRSkNCmxpYnJhcnkodGlkeXZlcnNlKQ0KbGlicmFyeShnZ3JpZGdlcykNCmxpYnJhcnkoZ2FtbHNzKQ0KbGlicmFyeShyc3RhbikNCmxpYnJhcnkoYnJtcykNCmBgYA0KDQohW10oTkJJQmF5ZXMxLnBuZykNCg0KKlTDoWMgZ2nhuqM6KiANCg0KKsSQaW5oIFRp4bq/biBUw6BpICgxLDIsMykqDQoqTMOqIMSQw7RuZyBOaOG6rXQgTmFtICgxLDMsNCkqDQoqVHLhuqduIFRydW5nIETFqW5nICgzKSoNCipMw6ogTmfhu41jIEto4bqjIE5oaSAoMSw0LDUpKg0KDQpHaGkgY2jDujogxJDDs25nIGfDs3AgY+G7p2EgY8OhYyB0aMOgbmggdmnDqm46DQoNCigxKSBCw6BpIGdp4bqjbmcgbMO9IHRodXnhur90DQooMikgTcO0IGjDrG5oIFBvaXNzb24gdsOgIE5CSSBj4buVIMSRaeG7g24NCigzKSBTb+G6oW4gdsOgIGhvw6BuIHRoaeG7h24gU1RBTiBjb2RlcyBjaG8gbcO0IGjDrG5oIE5CSSBCYXllcw0KKDQpIEtoYWkgdGjDoWMgbcO0IGjDrG5oIHbDoCBTdXkgZGnhu4VuIEJheWVzDQooNSkgxJDhu5MgaOG7jWEgdGjhu5FuZyBrw6oNCg0KIyBHaeG7m2kgdGhp4buHdQ0KDQpUaMOibiBjaMOgbyBjw6FjIGLhuqFuLCDEkcOieSBsw6AgYsOgaSB0aOG7sWMgaMOgbmggdGjhu6kgNyB0cm9uZyBk4buxIMOhbiBCYXllcyBmb3IgVmlldG5hbSB24bubaSBz4buxIGfDs3AgbeG6t3QgY+G7p2EgdG/DoG4gdGjhu4MgY29yZXRlYW0gY+G7p2EgcHJvamVjdC4gTeG7pWMgdGnDqnUgY+G7p2EgbmjDs20gQkFWIGzDoCBwaOG7lSBj4bqtcCB24buBIHBoxrDGoW5nIHBow6FwIHRo4buRbmcga8OqIHRoZW8gdHLGsOG7nW5nIHBow6FpIEJheWVzIG5o4bqxbSB0aGF5IHRo4bq/IGhvw6BuIHRvw6BuIG5o4buvbmcgY8O0bmcgY+G7pSB0cnV54buBbiB0aOG7kW5nLiDEkOG7kWkgdMaw4bujbmcgY+G7p2EgY2jDum5nIHTDtGkgbMOgIGPDoWMgYuG6oW4gYsOhYyBzxKkgdsOgIHNpbmggdmnDqm4geSBraG9hLiANCg0KxJDDonkgY8WpbmcgbMOgIG3hu5l0IGLDoGkgZ2nhuqNuZyB24bubaSB0aGFtIHbhu41uZyBjYW8gbmjhuqV0LHbDrCBjaMO6bmcgdMO0aSBz4bq9IMSR4buBIGPhuq1wIMSR4buTbmcgdGjhu51pIG5oaeG7gXUgduG6pW4gxJHhu4EgdMawxqFuZyDEkeG7kWkgcGjhu6ljIHThuqFwLCBiYW8gZ+G7k206IFBow6JuIHTDrWNoIGJp4bq/biBz4buRIMSR4bq/bSAoY291bnQgZGF0YSkgYuG6sW5nIHBow6JuIHBo4buRaSBOaOG7iyB0aOG7qWMgw6JtLCBtaXhlZCBtb2RlbCBjw7MgY2jhu6lhIHJhbmRvbSBlZmZlY3QsIHbDoCBI4buTaSBxdXkgdHV54bq/biB0w61uaCDEkWEgYmnhur9uIChuaOG7r25nIGLDoGkgZ2nhuqNuZyB0csaw4bubYyBraWEgY2jhu4kgZOG7q25nIGzhuqFpIOG7nyBwaMOibiB0w61jaCDEkcahbiBiaeG6v24pLg0KDQpDaMO6bmcgdMO0aSBnaeG6oyDEkeG7i25oIHLhurFuZyBjw6FjIGLhuqFuIMSRw6MgYuG6r3QgxJHhuqd1IHF1ZW4gduG7m2kgY+G6pXUgdHLDumMgbmfDtG4gbmfhu68gU1RBTiwgcXV5IHRyw6xuaCBjaHV54buDbiBnaeG6oyB0aHV54bq/dCBuZ2hpw6puIGPhu6l1IHRow6BuaCBtw7QgaMOsbmggQmF5ZXMsIGtoYWkgdGjDoWMgcGjDom4gcGjhu5EgaOG6rXUgxJHhu4tuaCAuTmjhu69uZyBwaOG6p24gbsOgeSBz4bq9IMSRxrDhu6NjIGdp4bqjbiBsxrDhu6NjDQoNCiMgVHLGsOG7nW5nIGjhu6NwIG1pbmggaOG7jWENCg0KYGBge3IsbWVzc2FnZSA9IEZBTFNFLHdhcm5pbmc9RkFMU0V9DQpsaWJyYXJ5KHRpZHl2ZXJzZSkNCmxpYnJhcnkoZ2dyaWRnZXMpDQoNCm15X3RoZW1lIDwtIGZ1bmN0aW9uKGJhc2Vfc2l6ZSA9IE5VTEwsIGJhc2VfZmFtaWx5ID0gInNhbnMiKSB7DQogIHRoZW1lX2J3KGJhc2Vfc2l6ZSA9IGJhc2Vfc2l6ZSwgYmFzZV9mYW1pbHkgPSBiYXNlX2ZhbWlseSkgKw0KICAgIHRoZW1lKA0KICAgICAgcGFuZWwuZ3JpZC5tYWpvciA9IGVsZW1lbnRfbGluZShjb2xvciA9ICJncmF5IiksDQogICAgICBwYW5lbC5ncmlkLm1pbm9yID0gZWxlbWVudF9ibGFuaygpLA0KICAgICAgcGFuZWwuYmFja2dyb3VuZCA9IGVsZW1lbnRfcmVjdChmaWxsID0gIiNmZmY2OGYiICksDQogICAgICBzdHJpcC5iYWNrZ3JvdW5kID0gZWxlbWVudF9yZWN0KGZpbGwgPSAiIzUxMDE2NiIsIGNvbG9yID0gTkEsIHNpemUgPTAuNSksDQogICAgICBzdHJpcC50ZXh0ID0gZWxlbWVudF90ZXh0KGZhY2UgPSAiYm9sZCIsIHNpemUgPSA4LCBjb2xvciA9ICJ3aGl0ZSIpLA0KICAgICAgbGVnZW5kLnBvc2l0aW9uID0gImJvdHRvbSIsDQogICAgICBsZWdlbmQuanVzdGlmaWNhdGlvbiA9ICJjZW50ZXIiLA0KICAgICAgbGVnZW5kLmJhY2tncm91bmQgPSBlbGVtZW50X2JsYW5rKCksDQogICAgICBsZWdlbmQubWFyZ2luID0gbWFyZ2luKDAuNSwwLjUsMC41LDAuNSkNCiAgICApDQp9DQpgYGANCg0KQuG7mSBk4buvIGxp4buHdSBuw6B5IMSRxrDhu6NjIHRodSB0aOG6rXAgdOG7qyBt4buZdCBuZ2hpw6puIGPhu6l1IHRyw6puIDU5IGLhu4duaCBuaMOibiBt4bqvYyBi4buHbmggxJHhu5luZyBraW5oLiBDw6FjIGLhu4duaCBuaMOibiDEkcaw4bujYyB0aGVvIGTDtWkgc+G7kSBs4bqnbiBi4buLIMSR4buZbmcga2luaCB0cm9uZyB2w7JuZyA4IHR14bqnbiAoYmFzZSkuIFNhdSDEkcOzIGLhu4duaCBuaMOibiDEkcaw4bujYyBwaMOibiBjaGlhIG5n4bqrdSBuaGnDqm4gdsOgbyBuaMOzbSDEkWnhu4F1IHRy4buLICh0cmVhdG1lbnQtIGTDuW5nIHRodeG7kWMgUHJvZ2FiaWRlICkgdsOgIG5ow7NtIGdp4bqjIGTGsOG7o2MgKHBsYWNlYm8pOyBT4buRIGzhuqduIGLhu4sgxJHhu5luZyBraW5oIMSRxrDhu6NjIGdoaSBuaOG6rW4gbeG7l2kgMiB0deG6p24gdHJvbmcgIDggdHXhuqduIHRp4bq/cCB0aGVvLiBOaMawIHbhuq15IG3hu5dpIGLhu4duaCBuaMOibiBz4bq9IMSRxrDhu6NjIGto4bqjbyBzw6F0IOG7nyA0IHRo4budaSDEkWnhu4NtIChwZXJpb2QpIHbDoCBk4bqhbmcgZOG7ryBsaeG7h3UgbsOgeSDEkcaw4bujYyBn4buNaSBsw6AgZOG7ryBsaeG7h3UgdMOhaSDEkW8gbMaw4budbmcgKHJlcGVhdGVkIG1lYXN1cmVtZW50IGRhdGEpLiBEYXRhc2V0IMSRw6MgxJHGsOG7o2Mgc+G7rSBk4bulbmcgYuG7n2kgbmhp4buBdSB0w6FjIGdp4bqjIG5oxrAgVGhhbGwgLFZhaWwgWzE5OTBdLCBCcmVzbG93IHbDoCBDbGF5dG9uWzE5OTNdLCBMZWUgdsOgIE5lbGRlciBbMTk5NiwgMjAwMF0gxJHhu4MgbWluaCBo4buNYSBjaG8gbcO0IGjDrG5oIGjhu5NpIHF1eSBo4buXbiBo4bujcCB2w6AgY291bnQgZGF0YS4gDQoNCsSQw6J5IGzDoCBt4buZdCB0aOG7rSBuZ2hp4buHbSBsw6JtIHPDoG5nIHRpw6p1IGJp4buDdSwgdHJvbmcgxJHDsyBjw6J1IGjhu49pIG5naGnDqm4gY+G7qXUgbMOgIGto4bqjbyBzw6F0IHTDoWMgxJHhu5luZyBj4bunYSBt4buZdCBsb+G6oWkgdGh14buRYyBYIGzDoG0gdGhheSDEkeG7lWkgKMSRxrDhu6NjIGvDrCB24buNbmcgbMOgIHTDrWNoIGPhu7FjKSBt4buZdCBr4bq/dCBj4bulYyB44bqldSBj4bunYSBt4buZdCBi4buHbmggbMO9LiBUaGnhur90IGvhur8gbmdoacOqbiBj4bupdSB0aMO0bmcgZOG7pW5nIGzDoCBwaMOibiBjaGlhIG5n4bqrdSBuaGnDqm4gdsOgIMSR4buTbmcgbmjhuqV0IG3hu5l0IG3huqt1IG4gYuG7h25oIG5ow6JuIHbDoG8gMiBwaMOibiBuaMOzbSAobmjDoW5oKSB0cuG7iyBsaeG7h3UsIG3hu5l0IGNobyBkw7luZyB0aHXhu5FjIFgsIG3hu5l0IHPhu60gZOG7pW5nIGdp4bqjIGTGsOG7o2MgKHBsYWNlYm8pLiBL4bq/dCBj4bulYyAob3V0Y29tZSkgWSBz4bq9IMSRxrDhu6NjIHRoZW8gZMO1aSBrw6lvIGTDoGkgcXVhIG5oaeG7gXUgdGjhu51pIMSRaeG7g20gdHJvbmcgc3Xhu5F0IHF1w6EgdHLDrG5oIMSRaeG7gXUgdHLhu4suIE3hu5l0IG3DtCBow6xuaCBo4buTaSBxdXkgxJFhIGJp4bq/biB0cm9uZyDEkcOzIHBow6JuIG5ow7NtIMSRaeG7gXUgdHLhu4sgKFggdnMgcGxhY2VibykgbMOgIGJp4bq/biBz4buRIGNow61uaCwga8OobSB0aGVvIG5o4buvbmcgaGnhu4dwIGJp4bq/biBuaMawIHRo4budaSBnaWFuLCB0deG7lWksIHTDrG5oIHRy4bqhbmcgY8ahIGLhuqNuIChiYXNlbGluZSnigKYgIHbDoCB54bq/dSB04buRIG5n4bqrdSBuaGnDqm4gY8OhIHRo4buDL3Ro4bupIGLhuq1jIGtow6FjLiBTdXkgZGnhu4VuIHRo4buRbmcga8OqIMSRxrDhu6NjIHRo4buxYyBoaeG7h24gdHLDqm4gdGhhbSBz4buRIGjhu5NpIHF1eSBjaG8gWC4NCg0KVHJvbmcgdHLGsOG7nW5nIGjhu6NwIG7DoHksIGvhur90IGPhu6VjIGzDom0gc8OgbmcgWSBsw6Agc+G7kSBjxqFuIMSR4buZbmcga2luaCDEkcaw4bujYyBnaGkgbmjhuq1uIHTDrW5oIMSR4bq/biBt4buXaSBt4buRYyB0aOG7nWkgZ2lhbiAoaz00KSwgeeG6v3UgdOG7kSBjYW4gdGhp4buHcCBjw7MgMiBsZXZlbHMgbMOgIHRodeG7kWMgcHJvZ2FiaWRlIHNvIHbhu5tpIHBsYWNlYm8uIEdp4bqjIHRodXnhur90IG5naGnDqm4gY+G7qXUgbMOgIHRodeG7kWMgUHJvZ2FiaWRlIGPhuqNpIHRoaeG7h24gY8OzIMO9IG5naMSpYSAobMOgbSBnaeG6o20gc+G7kSBjxqFuIMSR4buZbmcga2luaCkgc28gduG7m2kgUGxhY2Viby4gQsOgaSB0b8OhbiBj4bunYSBjaMO6bmcgdGEgY8OybiB4w6l0IHRow6ptIDIgaGnhu4dwIGJp4bq/biBz4buRIHRp4buBbSDhuqluIGtow6FjIMSRw7MgbMOgIDogVMOsbmggdHLhuqFuZyBjxqFuIMSR4buZbmcga2luaCAobGJhc2UgOiBsb2dhcml0IGPhu6dhIHThuqduIHN14bqldCDEkeG7mW5nIGtpbmggdHLGsOG7m2Mga2hpIHBow6JuIG5ow7NtICkgdsOgIFRo4budaSBnaWFuICh2aXNpdCA6IGzGsHUgw70sIHRyb25nIHRow60gbmdoaeG7h20gZ+G7kWMsIHRo4budaSBnaWFuIMSRxrDhu6NjIMSRbyBuaMawIG3hu5l0IGJp4bq/biBsacOqbiB04bulYywgdMOtbmggYuG6sW5nIHRow6FuZywgdsOgIGPDsyA0IHRo4budaSDEkWnhu4NtIGto4bqjbyBzw6F0LCB0dXkgbmhpw6puIGNow7puZyB0w7RpIGhvw6FuIGNodXnhu4NuIGJp4bq/biBz4buRIG7DoHkgdGjDoG5oIG5o4buvbmcgdHLhu41uZyBz4buRIHTGsMahbmcgcGjhuqNuIDogY29udHJhc3Qgd2VpZ2h0cykgxJHhu4MgxJHGoW4gZ2nhuqNuIGjDs2EgbcO0IGjDrG5oLiBWaeG7h2MgeMOpdCB0aMOqbSAyIGJp4bq/biBz4buRIG7DoHkgbmjhurFtIMSR4buDIGhp4buHdSBjaOG7iW5oIHNhaSBs4buHY2ggKG7hur91IGPDsykgZG8gxJHhu5kgbuG6t25nIGPDoSB0aOG7gyBsw6puIHN1eSBkaeG7hW4gdGjhu5FuZyBrw6ogduG7gSBoaeG7h3UgcXXhuqMgxJFp4buBdSB0cuG7iywgY8WpbmcgbmjGsCBoaeG7h3Ug4bupbmcgcGjhu6UgY+G7p2EgdGjhu51pIGdpYW4uDQoNCmBgYHtyLG1lc3NhZ2UgPSBGQUxTRSx3YXJuaW5nPUZBTFNFfQ0KZGF0PXJlYWQuY3N2KCJodHRwczovL3Jhdy5naXRodWIuY29tL3ZpbmNlbnRhcmVsYnVuZG9jay9SZGF0YXNldHMvbWFzdGVyL2Nzdi9NQVNTL2VwaWwuY3N2IikNCm5hbWVzKGRhdClbMl0gPC0gInNlaSINCmRhdCR2aXNpdCA8LSAoMipkYXQkcGVyaW9kLTUpLzEwDQoNCmRhdCRwZXJpb2Q9YXMuZmFjdG9yKGRhdCRwZXJpb2QpDQpkYXQkc3ViamVjdD1hcy5mYWN0b3IoZGF0JHN1YmplY3QpDQoNCmRhdCU+JWhlYWQoKSU+JSBrbml0cjo6IGthYmxlKCkNCmBgYA0KDQpOaMawbmcgdHLGsOG7m2Mga2hpIMSRaSB2w6BvIGLDoGkgdG/DoW4sIGNow7puZyB0w7RpIG114buRbiBiw6BuIG3hu5l0IGNow7p0IHbhu4EgbMO9IHRodXnhur90IHjDoWMgc3XhuqV0IHbDoCBtw7QgaMOsbmggY2hvIGJp4bq/biBz4buRIMSR4bq/bSAoY291bnQgZGF0YSkNCg0KVGjhu5FuZyBrw6ogbcO0IHThuqMgY2hvIHRo4bqleSBPdXRjb21lIChTZWl6dXJlIGNvdW50KSBjw7MgcGjDom4gcGjhu5FpIGtow7RuZyBiw6xuaCB0aMaw4budbmcgOiBD4bulIHRo4buDLCBuw7MgY8OzIHRydW5nIGLDrG5oIG5o4buPIGjGoW4gU2QsIHNrZXduZXNzIHLhuqV0IGNhby4NCg0KYGBge3J9DQpwc3ljaDo6ZGVzY3JpYmUoZGF0JHNlaSkNCg0Kc2lnbWE9c2QoZGF0JHNlaSkNCm11PW1lYW4oZGF0JHNlaSkNCg0KcmVhbHA9ZGF0JT4lZ2dwbG90KGFlcyh4PXNlaSkpKw0KICBnZW9tX2RlbnNpdHkoYWxwaGE9LjgsZmlsbD0iI2M1MDdmZiIpKw0KICBnZW9tX2hpc3RvZ3JhbShhZXMoeT0uLmRlbnNpdHkuLixmaWxsPS4uZGVuc2l0eS4uKSxjb2xvdXI9ImJsYWNrIixhbHBoYT0wLjcsc2hvdy5sZWdlbmQgPSBGLGJpbndpZHRoID0gMSkrDQogIHRoZW1lX2J3KCkrc2NhbGVfZmlsbF9ncmFkaWVudChsb3c9IiNmZjljMDciLGhpZ2g9IiNmZjA3MmYiKSsNCiAgZ2d0aXRsZSgiUmVhbCBkYXRhIikNCg0KcmVhbHANCmBgYA0KDQpOZ3V5w6puIG5ow6JuIHbDrCDEkcOieSBsw6AgbeG7mXQgc+G7kSDEkeG6v20gKGNvdW50IGRhdGEpLg0KDQpUcm9uZyB5IGjhu41jIGzDom0gc8OgbmcsIG5oaeG7gXUgxJHhuqFpIGzGsOG7o25nIMSRxrDhu6NjIGto4bqjbyBzw6F0IG5oxrAgbeG7mXQgc+G7kSDEkeG6v20gLyBiaeG6v24gc+G7kSBy4budaSBy4bqhYyAoY291bnQgZGF0YSwgZGlzY3JldGUgdmFyaWFibGVzKS4gQ2jDum5nIHRoxrDhu51uZyBtYW5nIMO9IG5naMSpYSB04bqnbiBzdeG6pXQgKHPhu5EgbOG6p24pIG3hu5l0IGJp4bq/biBj4buRL2hp4buHbiB0xrDhu6NuZyDEkcaw4bujYyBnaGkgbmjhuq1uIHRyb25nIG3hu5l0IGtob+G6o25nIGtow7RuZyBnaWFuIHbDoCB0aOG7nWkgZ2lhbiB4w6FjIMSR4buLbmgsIFRow60gZOG7pSA6DQoNCisgSGnhu4duIHTGsOG7o25nIHNpbmggbMO9IDogbmjhu4twIHRpbSwgbmjhu4twIGjDtCBo4bqlcCwgbmh1IMSR4buZbmcgcnXhu5l0LCBz4buRIGPGoW4gZ8OyIHThu60gY3VuZywgDQoNCisgQmnhur9uIGPhu5EgYuG7h25oIGzDvSA6IFPhu5EgY8ahbiBoZW4ga+G7i2NoIHBow6F0LCBz4buRIGzhuqduIG5nxrBuZyB0aOG7nyBraGkgbmfhu6csIHPhu5EgY8ahbiDEkeG7mW5nIGtpbmgsIHPhu5EgdHLGsOG7nW5nIGjhu6NwIHThu60gdm9uZywgc+G7kSBs4bqnbiBuaOG6rXAgdmnhu4duIGPhuqVwIGPhu6l14oCmDQoNCisgS+G6v3QgcXXhuqMgY+G7p2EgcGjDqXAgxJHhur9tIDogU+G7kSBsxrDhu6NuZyB04bq/IGLDoG8gYuG6oWNoIGPhuqd1LCBz4buRIGvDrSBzaW5oIHRyw7luZywgc+G7kSBjb24sIHPhu5EgxJHGoW4gduG7iyBjxqEgcXVhbiBi4buLIHThu5VuIHRoxrDGoW5n4oCmIA0KDQorIEjDoG5oIHZpIDogc+G7kSBs4bqnbiB0aMSDbSBraMOhbSBiw6FjIHPEqSwgbeG7qWMgdGnDqnUgdGjhu6UgdGh14buRYyBsw6EvYmlhIHLGsOG7o3UsIHPhu5EgbOG6p24gZMO5bmcgdGh14buRYyBj4bqvdCBjxqFuIGhlbg0KDQorIEhv4bq3YyBt4buZdCB0aGFuZyDEkW8gcuG7nWkgcuG6oWMgOiDEkWnhu4NtIEdsYXNnb3csIG3hu6ljIMSR4buZIMSRYXUsIA0KDQorIFRo4budaSBnaWFuIGPFqW5nIGPDsyB0aOG7gyDEkcaw4bujYyB4w6l0IG5oxrAgYmnhur9uIHPhu5EgxJHhur9tLCBraGkgY2jDum5nIGzDoCBz4buRIG5ndXnDqm4sIHRow60gZOG7pTogU+G7kSBuZ8OgeSBu4bqxbSB2aeG7h24sIMSR4buZIGTDoGkgZ2nhuqVjIG5n4bunICwuLi4gDQoNCsSQ4bq3YyDEkWnhu4NtIGPhuqduIGzGsHUgw70gduG7m2kgY291bnQgZGF0YSBsw6AgcGjDom4gYuG7kSBk4buvIGxp4buHdSB0aMaw4budbmcgaGF5IGLhu4sgbOG7h2NoIHBo4bqjaSwgcGjGsMahbmcgc2FpIHTEg25nIGtoaSBnacOhIHRy4buLIHRydW5nIGLDrG5oIHTEg25nLiAgDQoNCiMgQmnhur9uIHPhu5EgxJHhur9tOiB04burIMSR4bqhaSBsxrDhu6NuZyDEkeG6v24gYmnhur9uIG5n4bqrdSBuaGnDqm4NCg0KVGnhur9wIHRoZW8sIGNow7puZyB0w7RpIG114buRbiBwaMOibiBiaeG7h3QgZ2nhu69hIMSR4bqhaSBsxrDhu6NuZyB2w6AgYmnhur9uIHPhu5EgdHJvbmcgbmdoacOqbiBj4bupdSwgdsOsIGPDsyBz4buxIGtow6FjIG5oYXUgZ2nhu69hIDIga2jDoWkgbmnhu4dtIG7DoHkNCg0KIVtdKG1lYXN1cmUucG5nKQ0KDQpYw6FjIMSR4buLbmggxJHGsOG7o2MgYuG6o24gY2jhuqV0LCDDvSBuZ2jEqWEgdsOgIGPDoWNoIMSRbyDEkeG6v20gbmjhu69uZyDEkeG6oWkgbMaw4bujbmcgbMOgIG3hu5l0IGLGsOG7m2MgcXVhbiB0cuG7jW5nIGtoaSDEkeG6t3QgZ2nhuqMgdGh1eeG6v3QgdsOgIHRoaeG6v3Qga+G6vyBuZ2hpw6puIGPhu6l1LCDEkcOieSBsw6AgdHLDoWNoIG5oaeG7h20gY+G7p2EgbmfGsOG7nWkgYsOhYyBzxKkuIE5oxrBuZyBixrDhu5tjIHRp4bq/cCB0aGVvLCBxdWFuIHRy4buNbmcga2jDtG5nIGvDqW0sIMSRw7MgbMOgIGNodXnhu4NuIHThu6sgxJHhuqFpIGzGsOG7o25nIHNhbmcgYmnhur9uIHPhu5Eg4oCTdsOgIMSRw6J5IGzDoCBjw7RuZyB2aeG7h2MgY+G7p2EgY2h1ecOqbiB2acOqbiB0aOG7kW5nIGvDqi4gDQoNCk3hu5l0IGPDoWNoIGzDvSB0xrDhu59uZywgbmfGsOG7nWkgYsOhYyBzxKkgY8OzIMSR4bqneSDEkeG7pyBj4bqjIGtp4bq/biB0aOG7qWMgdGjhu5FuZyBrw6ogdsOgIHkgaOG7jWMsIHPhur0gY8O5bmcgbMO6YyDEkeG6o20gbmjhuq1uIGPhuqMgMiB2YWkgdHLDsiDEkeG7gyB04buxIG3DrG5oIMSRaSB04burIHF1YW4gc8OhdCDEkeG6v24gZ2nhuqMgdGh1eeG6v3QsIHThu6sgZ2nhuqMgdGh1eeG6v3QgxJHhur9uIHRow60gbmdoaeG7h20sIMSRby/EkeG6v20gxJHhuqFpIGzGsOG7o25nIHbDoCB04burIMSR4bqhaSBsxrDhu6NuZyBuw6B5IMSR4bq/biBiaeG6v24gc+G7kSB0aOG7kW5nIGvDqi4gDQoNCkPDoWMgbcO0IGjDrG5oIHRo4buRbmcga8OqIGzDoCBuaOG7r25nIGNvbiByb2JvdCB2w7QgdHJpLCBjaMO6bmcgaG/DoG4gdG/DoG4ga2jDtG5nICJuaMOsbiIgdGjhuqV5IMSR4bqhaSBsxrDhu6NuZywga2jDtG5nIGJp4bq/dCDEkeG6v24gw70gbmdoxKlhIGzDom0gc8OgbmcuIENow7puZyBjaOG7iSBiaeG6v3QgdOG7m2kgYmnhur9uIHPhu5Egbmfhuqt1IG5oacOqbi4gQ2jDum5nIGtow7RuZyDEkW8gbMaw4budbmcsIG3DoCBjaOG7iSBkw7luZyBow6BtIHRvw6FuIGjhu41jIMSR4buDIMaw4bubYyBsxrDhu6NuZyBnacOhIHRy4buLIGPhu6dhIGJp4bq/biBz4buRIGPhuqduIHjDqXQuIE3hu5dpIGdpw6EgdHLhu4sgbmjGsCB24bqteSBs4bqhaSBjw7MgeMOhYyBzdeG6pXQgaGnhu4duIGRp4buHbiBraMOhYyBuaGF1LiBOaMOgIHRo4buRbmcga8OqIGtow7RuZyBjw7Mga2nhur9uIHRo4bupYyB5IGjhu41jIGNodXnDqm4gbcO0biBz4bq9IGtow7RuZyBxdWFuIHTDom0gxJHhur9uIMO9IG5naMSpYSBzaW5oIGzDvSwgYuG7h25oIGjhu41jIG7hu69hLCBuaMawbmcgY2jhu4kgY+G7kSBn4bqvbmcgbcO0IHThuqMgxJHhurdjIHTDrW5oIHBow6JuIHBo4buRaSBj4bunYSBiaeG6v24gbmfhuqt1IG5oacOqbiBi4bqxbmcgcXV5IGx14bqtdCBwaMO5IGjhu6NwIChow6BtIG3huq10IMSR4buZIHjDoWMgc3XhuqV0KS5Dw6FjIGjDoG0gbsOgeSDEkcaw4bujYyDEkeG7i25oIG5naMSpYSBi4bqxbmcgY8OhYyB0aGFtIHPhu5EgdsOgIGdp4bqjIMSR4buLbmguIFThu6sgcGjDom4gcGjhu5FpLCBo4buNIGPDsyB0aOG7gyDGsOG7m2MgbMaw4bujbmcgdGhhbSBz4buRL2dpw6EgdHLhu4sgYuG6sW5nIG3DtCBow6xuaCwgdsOgIHN1eSBkaeG7hW4gdGjhu5FuZyBrw6ouIA0KDQpU4bubaSDEkcOieSwgY8OhYyBi4bqhbiBzaW5oIHZpw6puIHbDoCBiw6FjIHPEqSBjxaluZyBjw7MgdGjhu4Mgbmjhuq1uIHJhIHbhuqVuIMSR4buBIDogTmjhu69uZyBiaeG6v24gc+G7kSDEkeG6v20sIHLhu51pIHLhuqFjIGtow7RuZyB0aOG7gyDEkcaw4bujYyBtw7QgdOG6oyBi4bqxbmcgY8O5bmcgbeG7mXQgcGjDom4gcGjhu5FpIEdhdXNzaWFuIChjaHXhuqluLCBiw6xuaCB0aMaw4budbmcgKSBuaMawIG5o4buvbmcgYmnhur9uIHPhu5EgbGnDqm4gdOG7pWMuIE5oxrAgduG6rXkgY8OhYyBwaMawxqFuZyBwaMOhcCB44butIGzDvSBz4buRIGxp4buHdSBjw7MgZ2nhuqMgxJHhu4tuaCBwaMOibiBwaOG7kWkgY2h14bqpbiAodC10ZXN0LCBBTk9WQSwgT0xTIHJlZ3Jlc3Npb24sIGV0Yykga2hpIMOhcCBk4bulbmcgduG7m2kgZOG7ryBsaeG7h3UgY291bnQgZGF0YSBz4bq9IGtow7RuZyBjw7JuIHBow7kgaOG7o3AgbuG7r2EgKHRy4burIGtoaSBk4buvIGxp4buHdSDEkcOjIMSRxrDhu6NjIGNodXnhu4NuIMSR4buVaSkuIEtoaSBi4bqhbiDEkcOjIGdp4bqjIMSR4buLbmggc2FpIHbhu4EgcXV5IGx14bqtdCBwaMOibiBwaOG7kWksIG3DtCBow6xuaCBj4bunYSBi4bqhbiDEkcOjIGPDsyBiaWFzIHbDoCBuaGnhu4F1IG5ndXkgY8ahIGzDoCBzdXkgZGnhu4VuIHRo4buRbmcga8OqIGThu7FhIHRyw6puIG3DtCBow6xuaCBuw6B5IGPFqW5nIHNhaSBs4bqnbSBu4buRdC4NCg0KxJBp4buBdSDEkcOhbmcgdGnhur9jIMSRw7MgbMOgIGNoxrDGoW5nIHRyw6xuaCBUaOG7kW5nIGvDqiDhu58gxJHhuqFpIGjhu41jIFkgY2jhu4kgcXVhbiB0w6JtIMSR4bq/biBuaOG7r25nIGJp4bq/biBsacOqbiB04bulYywgbmjGsG5nIGNoxrBhIHRyYW5nIGLhu4sgY2hvIGPDoWMgc2luaCB2acOqbiBZIGtob2EgxJHhu6cga2nhur9uIHRo4bupYyDEkeG7gyB44butIHRyw60gbmjhu69uZyBiaeG6v24gcuG7nWkgcuG6oWMgdsOgIGJp4bq/biBz4buRIMSR4bq/bS4NCg0KTeG7mXQgbcOidSB0aHXhuqtuIGtow6FjLCDEkcOzIGzDoCBraGkgdGhpIG7hu5lpIHRyw7osIGNhbyBo4buNYywgbmfGsOG7nWkgdGEgY8OzIGThuqF5IGzDvSB0aHV54bq/dCB4w6FjIHN14bqldCB24buBIG5o4buvbmcgcXV5IGx14bqtdCBwaMOibiBwaOG7kWkgbmjGsCBCaW5vbWlhbCAoTmjhu4sgdGjhu6ljICksIEJlcm5vdWxsaSwgUG9pc3NvbiwgbmjGsG5nIGPDoWNoIGThuqF5IHRodeG6p24gdMO6eSBsw70gdGh1eeG6v3QgdsOgIHRvw6FuIGjhu41jIGtow7RuZyBjaG8gcGjDqXAgduG6rW4gZOG7pW5nIG5o4buvbmcgcXV5IGx14bqtdCBwaMOibiBwaOG7kWkgbsOgeSB2w6BvIHRo4buxYyB0aeG7hW4gbMOibSBzw6BuZyB2w6AgbmdoacOqbiBj4bupdSDigJMgc2luaCB2acOqbiBraMO0bmcgY8OzIGto4bqjIG7Eg25nIGxpw6puIGjhu4cgc+G7sSB24bqtdCwgaGnhu4duIHTGsOG7o25nIHbhu5tpIHF1eSBsdeG6rXQgcGjDom4gcGjhu5FpLCB2w6Aga2jDtG5nIGPDsyBraOG6oyBuxINuZyBjaHV54buDbiB04burIMSR4bqhaSBsxrDhu6NuZyB0aMOgbmggYmnhur9uIHPhu5Egbmfhuqt1IG5oacOqbi4gDQoNCk3hu5l0IHNhaSBs4bqnbSBjaOG6v3QgbmfGsOG7nWkga2jDoWMsIMSRw7MgbMOgIHZp4buHYyBk4bqheSBjaG8gc2luaCB2acOqbiB0aeG6v3AgY+G6rW4gduG6pW4gxJHhu4EgYuG6sW5nIHTGsCBkdXkgc28gc8OhbmgsIHbDoCBjw6FjIGPDtG5nIGPhu6UgbOG6oWMgaOG6rXUgbmjGsCBi4bqjbmcgY2jDqW8sIGtp4buDbSDEkeG7i25oIHBoaSB0aGFtIHPhu5HigKYgbmjGsG5nIGtow7RuZyBk4bqheSB24buBIG3DtCBow6xuaCBo4buTaSBxdXksIG5o4bqldCBsw6AgbcO0IGjDrG5oIHR1eeG6v24gdMOtbmggdOG7lW5nIHF1w6F0IChHTE0pLiANCg0KVuG7m2kgMiBs4buXIGjhu5VuZyBraeG6v24gdGjhu6ljIGPhuqMgduG7gSBYw6FjIHN14bqldCDhu6luZyBk4bulbmcsIE3DtCBow6xuaCBHTE0sIHNpbmggdmnDqm4gZ+G6p24gbmjGsCBob8OgbiB0b8OgbiBi4bqldCBs4buxYyB0csaw4bubYyBuaOG7r25nIGdp4bqjIHRodXnhur90IG5naGnDqm4gY+G7qXUgduG7m2kgYmnhur9uIHPhu5EgxJHhur9tLiBTaW5oIHZpw6puIGNo4buJIGPDsyAyIGzhu7FhIGNo4buNbiwgbMOgbSBsaeG7gXUgdsOgIHNhaSAoZMO5bmcga2nhu4NtIMSR4buLbmggdCBzbyBzw6FuaCwgduG7kW4gZOG7sWEgdHLDqm4gZ2nhuqMgxJHhu4tuaCBwaMOibiBwaOG7kWkgY2h14bqpbiksIGhv4bq3YyB0w6xtIG7GoWkgdHLDuiDhuqluIChjw6FjIHRlc3QgcGhpIHRoYW0gc+G7kSkuDQoNClRyb25nIGThu7Egw6FuIEJheWVzIGFuYWx5c2lzIGZvciBWaWV0bmFtLCBuaMOzbSBjaMO6bmcgdMO0aSBuaOG6r20gxJHhur9uIGPhuqMgMiBt4bulYyB0acOqdSA6IE3hu6VjIHRpw6p1IGNow61uaCBsw6Aga2h1eeG6v24ga2jDrWNoIGPDoWMgYuG6oW4gxJFpIHRoZW8gdHLGsOG7nW5nIHBow6FpIEJheWVzLCBuaMawbmcgbeG7mXQgbeG7pWMgdGnDqnUga2jDoWMgcXVhbiB0cuG7jW5nIGtow7RuZyBrw6ltLCDEkcOzIGzDoCBiw7kgxJHhuq9wIGNobyBjw6FjIGLhuqFuIMSR4buTbmcgbmdoaeG7h3Agbmjhu69uZyB0aGnhur91IGjhu6V0IHbhu4Ega2nhur9uIHRo4bupYyB2w6Aga+G7uSBuxINuZyB0aOG7kW5nIGvDqiDhu6luZyBk4bulbmcuIE3hu5dpIGLDoGkgdGjhu7FjIGjDoG5oIGzDoCBt4buZdCBjYXNlIHN0dWR5IGhvw6BuIGNo4buJbmggZOG6q24gZOG6r3QgY8OhYyBi4bqhbiB04burIGPDonUgaOG7j2ksIGdp4bqjIHRodXnhur90IMSR4bq/biBiaeG6v24gc+G7kSwgbcO0IGjDrG5oIHbDoCBzdXkgZGnhu4VuIGvhur90IHF14bqjLg0KDQojIFThu6sgcGjDom4gcGjhu5FpIG5o4buLIHRo4bupYyDEkeG6v24gTmjhu4sgdGjhu6ljIMOibQ0KDQpQaMOibiBwaOG7kWkgbmjhu4sgdGjhu6ljIHbDoCBCZXJub3VsbGkgbMOgIGfhu5FjIHLhu4UgY+G7p2EgbeG7jWkgcXV5IGx14bqtdCBwaMOibiBwaOG7kWkgZMOgbmggY2hvIGJp4bq/biBy4budaSBy4bqhYy4gDQoNClBow6JuIHBo4buRaSBuaOG7iyB0aOG7qWMgbcO0IHThuqMgYmnhur9uIHPhu5Egbmfhuqt1IG5oacOqbiB4IG5oxrAgU+G7kSBs4bqnbiB0aMOgbmggY8O0bmcgdHLDqm4gbiBs4bqnbiB0aOG7rSBuZ2hp4buHbSBs4bq3cCBs4bqhaSDEkeG7mWMgbOG6rXAgKG4geMOhYyDEkeG7i25oKSwgYmnhur90IHLhurFuZyB4w6FjIHN14bqldCB0aMOgbmggY8O0bmcgbMOgIHA6DQoNCiQkUCh4fG4scCk9XGJpbm9te259e3h9cF57eH0oMS1wKV57bi14fSQkDQoNClBow6JuIHBo4buRaSBuaOG7iyB0aOG7qWMgw6JtIDoNCg0KROG6oW5nIGNow61uaCB04bqvYyBj4bunYSBwaMOibiBwaOG7kWkgbmjhu4sgdGjhu6ljIMOibSBtw7QgdOG6oyBz4buRIGzhuqduIHRo4butIG5naGnhu4dtIEJlcm5vdWxsaSDEkeG7mWMgbOG6rXAsIGzhurdwIGzhuqFpIChiaeG6v24gbmfhuqt1IG5oacOqbiB4KSDEkeG7gyBjw7MgdGjhu4MgxJHhuqF0IMSRxrDhu6NjIG3hu5l0IHPhu5EgbOG6p24gdGjDoG5oIGPDtG5nIHjDoWMgxJHhu4tuaCAobOG6p24gdGjDoG5oIGPDtG5nIHRo4bupIHIpLCBiaeG6v3QgcuG6sW5nIHjDoWMgc3XhuqV0IHRow6BuaCBjw7RuZyBsw6AgcA0KDQokJFAoeHxwLHIpPSBcYmlub217eC0xfXtyLTF9cF57cn0oMS1wKV57eC1yfSQkDQpW4bubaSB4IHbDoCByIMSR4buBdSBsw6Agc+G7kSBuZ3V5w6puLCDEkWnhu4F1IGtp4buHbiB4ID49IHIsIHbDoCBwIGzDoCBt4buZdCB4w6FjIHN14bqldCB0cm9uZyBraG/huqNuZyAoMCwxKS4NCg0KROG6oW5nIHRo4bupIGhhaSBj4bunYSBwaMOibiBwaOG7kWkgbmjhu4sgdGjhu6ljIMOibSBsw6AgbeG7mXQgYmnhur9uIHRo4buDLCB0cm9uZyDEkcOzIGJp4bq/biBuZ+G6q3Ugbmhpw6puIHggbMOgIHPhu5EgbOG6p24gdGjhu60gbmdoaeG7h20gdGjhuqV0IGLhuqFpIHRyxrDhu5tjIGtoaSBs4bqnbiB0aMOgbmggY8O0bmcgdGjhu6kgciB44bqjeSByYS4gDQokJFAoeHxwLHIpPSBcYmlub217eCtyLTEpfXt4fXBee3J9KDEtcClee3h9JCQNClbhu5tpIHggbMOgIHPhu5Egbmd1ecOqbiBkxrDGoW5nIGhv4bq3Yz0wIHbDoCBraMO0bmcgY8OzIHLDoG5nIGJ14buZYyB24buBIMSRaeG7gXUga2nhu4duIGdp4buvYSByIHbDoCB4LiBOaMawIHbhuq15LCBjaMOtbmggZOG6oW5nIGJp4bq/biB0aOG7gyBuw6B5IG3hu5tpIHRo4buxYyBz4buxIGjhu691IGThu6VuZyDEkeG7gyBtw7QgdOG6oyBiaeG6v24gc+G7kSDEkeG6v20gdsOsIHggY8OzIHRo4buDIG5o4bqtbiBi4bqldCBj4bupIGdpw6EgdHLhu4sgc+G7kSBuZ3V5w6puIG7DoG8ga2jDtG5nIMOibSwgdsOgIHIgY8Wpbmcga2jDtG5nIGLhu4sgcsOgbmcgYnXhu5ljIMSRaeG7gXUga2nhu4duIG7DoG8gY+G6oy4NCg0KVMOqbiBn4buNaSAibmjhu4sgdGjhu6ljIMOibSBjw7MgduG6uyBiw60g4bqpbiA/IFRo4buxYyByYSBuw7MgcGjDoXQgeHXhuqV0IHThu6sgdHLGsOG7nW5nIGjhu6NwIHRhIMOhcCBk4bulbmcgcGjDom4gcGjhu5FpIG5o4buLIHRo4bupYyB24bubaSBz4buRIG3FqSDDom0gOg0KDQokJFxiaW5vbXt4K3ItMX17eH1cIHBee3J9KDEtcClee3h9ID0gXGJpbm9tey1yfXt4fVwgcF57cn0ocC0xKV57eH0kJA0KDQpUaGF5IHbDrCBiaeG7g3UgZGnhu4VuIHBow6JuIHBo4buRaSBuaOG7iyB0aOG7qWMgw6JtIGLhurFuZyAyIHRoYW0gc+G7kSByIHbDoCBwLCB0YSBjw7MgdGjhu4Mgc+G7rSBk4bulbmcga2jDoWkgbmnhu4dtIGtow6FjIGzDoCA6IHRoYW0gc+G7kSB24buLIHRyw60gdHJ1bmcgdMOibSAoTXUpIHbDoCBwaMawxqFuZyBzYWkgKHNpZ21hXjIpLCBxdWFuIGjhu4cgaG/DoW4gY2h1eeG7g24gZ2nhu69hIE11LCBTaWdtYV4yIHbDoCByLHAgbmjGsCBzYXUgOg0KDQoNCiQkXG11ID0gXGZyYWMge3IoMS1wKX17cH0kJA0KDQokJFZhciA9IFxzaWdtYV4yID0gXGZyYWN7cigxLXApfXtwXjJ9ID0gXG11ICsgXGZyYWN7MX17cn1cbXVeMiQkDQpE4bqhbmcgbsOgeSBjw7MgdMOqbiBn4buNaSBsw6AgTmVnYXRpdmUgYmlub21pYWwgdHlwZSAyIGhv4bq3YyBxdWFkcmF0aWMgbmVnYXRpdmUgYmlub21pYWwuIFPhu5EgMiBjaMOtbmggbMOgIGLhuq1jIG3FqSAyIGNobyBNdSB0cm9uZyBjw7RuZyB0aOG7qWMgdMOtbmggdmFyaWFuY2UuDQoNCiMgQmnhu4duIGx14bqtbiB24buBIHBow6JuIHBo4buRaSBjaG8gYsOgaSB0b8OhbiBtaW5oIGjhu41hDQoNClRy4bufIGzhuqFpIGLDoGkgdG/DoW4gY+G7p2EgY2jDum5nIHRhOg0KDQpOaOG7r25nIGzhuq1wIGx14bqtbiBzYXUgxJHDonkg4bunbmcgaOG7mSBjaG8gdmnhu4djIGNo4buNbiBo4buNIHBow6JuIHBo4buRaSBwaMO5IGjhu6NwIGNobyBtw7QgaMOsbmggOg0KDQpW4buBIG3hurd0IGzDom0gc8OgbmcgOiBPdXRjb21lIGPhu6dhIGNow7puZyB0YSBjw7MgYuG6o24gY2jhuqV0IGzDoCB04bqnbiBzdeG6pXQgcGjDoXQgc2luaCBt4buZdCBiaeG6v24gY+G7kSBsw6JtIHPDoG5nIChjxqFuIMSR4buZbmcga2luaCksIMSRxrDhu6NjIGto4bqjbyBzw6F0IHRyb25nIG3hu5l0IGtob+G6o25nIHRo4budaSBnaWFuIHjDoWMgxJHhu4tuaCwgdHJvbmcgbeG7mXQga2jDtG5nIGdpYW4geMOhYyDEkeG7i25oIChjxqEgdGjhu4MgbeG7l2kgYuG7h25oIG5ow6JuKS4gDQoNCk3hu5l0IHPhu5Egw70ga2jDoWMgY8OzIHRo4buDIMSRxrDhu6NjIHBow6F0IGJp4buDdSwgbmjGsDogTeG7l2kgY8ahbiDEkeG7mW5nIGtpbmggbMOgIMSR4buZYyBs4bqtcCAodGjhu7FjIHbhuq15ICEpLCBuZ+G6q3Ugbmhpw6puIChnaeG6oyDEkeG7i25oKSB2w6AgY8OzIHjDoWMgc3XhuqV0IGjhurFuZyDEkeG7i25oIHRoZW8gdGjhu51pIGdpYW4gKGdp4bqjIMSRaW5oKS4gTmjhu69uZyDEkeG6t2MgdMOtbmggbsOgeSB0aOG7j2EgxJFhIHPhu5EgxJFp4buBdSBraeG7h24gdsOgIHTDrW5oIGNo4bqldCBj4bunYSBt4buZdCBwaMOibiBwaOG7kWkgUG9pc3NvbiBkw6BuaCBjaG8gYmnhur9uIHPhu5EgxJHhur9tLCBuaMawbmcga2jDtG5nIHBow7kgaOG7o3AgduG7m2kgcGjDom4gcGjhu5FpIEdhdXNzaWFuLiBEbyDEkcOzLCBwaMOibiBwaOG7kWkgUG9pc3NvbiBsw6AgZ2nhuqMgxJHhu4tuaCBiYW4gxJHhuqd1Lg0KDQpLaGkgY8OhYyBi4bqhbiBtdeG7kW4geGVtIHjDqXQg4bqjbmggaMaw4bufbmcgY+G7p2EgY8OhYyBiaeG6v24gxJHhu5ljIGzhuq1wICB4MSx4MiwuLnhuIGzDqm4gYmnhur9uIHBo4bulIHRodeG7mWMgZOG6oW5nIHPhu5EgxJHhur9tIHkgdGjDrCBjaMO6bmcgdGEgdGjGsOG7nW5nIGhheSBuZ2jEqSBuZ2F5IMSR4bq/biBtw7QgaMOsbmggaOG7k2kgcXV5IFBvaXNzb24uIE3DtCBow6xuaCBQb2lzc29uIGPDsyBt4buZdCBnaeG6oyDEkeG7i25oIChhc3N1bXB0aW9uKSBy4bqldCBxdWFuIHRy4buNbmcgbMOgOiBnacOhIHRy4buLIG1lYW4gdsOgIHZhcmlhbmNlIGLhurFuZyBuaGF1LiBUdXkgbmhpw6puIGThu68gbGnhu4d1IHRo4buxYyB04bq/IG3DoCBjaMO6bmcgdGEgdGh1IHRo4bqtcCDDrXQga2hpIHRow7VhIG3Do24gZ2nhuqMgxJHhu4tuaCBuw6B5LiBLaGkgdmFyaWFuY2UgPiBtZWFuIGNow7puZyB0YSBn4buNaSDEkcOzIGzDoCBvdmVyLWRpc3BlcnNpb24gKHThuqFtIGThu4tjaDogcXXDoSBwaMOibiB0w6FuKSwgbmfGsOG7o2MgbOG6oWksIHZhcmlhbmNlIDwgbWVhbiB0aMOsIGfhu41pIGzDoCB1bmRlci1kaXNwZXJzaW9uLg0KDQpN4buZdCBiaeG6v24gbmfhuqt1IG5oacOqbiB5IHbhu5tpIGdpw6EgdHLhu4sga8OsIHbhu41uZyDOvCB0dcOibiB0aGVvIHBow6JuIGLhu5EgcG9pc3NvbiBz4bq9IGPDsyBow6BtIHByb2JhYmlsaXR5IG1hc3MgZnVuY3Rpb24gKGjDoG0ga2jhu5FpIHjDoWMgc3XhuqV0KSBuaMawIHNhdToNCg0KJCRQXGxlZnQgKCBZID0geSBccmlnaHQgKSA9IFxmcmFje2Veey1cbXUgfVxtdV55fXt5IX0kJA0KQ8OzIHRo4buDIG5o4bqtbiByYSBy4bqxbmcgcGjDom4gcGjhu5FpIHBvaXNzb24gY2jhu4kgY8OzIDEgdGhhbSBz4buRIGzDoCBtZWFuIChtdSkuIFRyb25nIGjhu5NpIHF1eSBwb2lzc29uLCBtdSDEkcaw4bujYyBtw7QgdOG6oyBi4bqxbmcgbeG7mXQgaMOgbSBz4buRIG3FqSBj4bunYSBjw6FjIGJp4bq/biDEkeG7mWMgbOG6rXAgdHV54bq/biB0w61uaCAobGluZWFyIHByZWRpY3RvcikuIExpbmsgZnVuY3Rpb24gxJHGsOG7o2Mgc+G7rSBk4bulbmcgdHJvbmcgaOG7k2kgcXV5IHBvaXNzb24gbMOgIGxvZyBsaW5rLiANCg0KJCRcbXUgPSBlXlxldGEkJA0KJCRcZXRhID0gXGJldGFfMCArIFxiZXRhXzEgeF8xICsgXGJldGFfMiB4XzIgKyAuLi4rIFxiZXRhX24geF9uJCQNCg0KxJDhuqd1IHRpw6puIGNow7puZyB0YSBz4bq9IGTDuW5nIG3DtCBow6xuaCBwb2lzc29uIGtow7RuZyBjw7MgcmFuZG9tIGVmZmVjdCDEkeG7gyBraeG7g20gdHJhIMSR4buZIHBow6JuIHTDoW4gY+G7p2EgZOG7ryBsaeG7h3UgKGRpc3BlcnNpb24pLg0KDQpgYGB7cn0NCmdsbS5wb2lzIDwtICBnbG0oc2VpIH4gbGJhc2UgKyB0cnQgKyB2aXNpdCAsIGZhbWlseSA9ICJwb2lzc29uIiwgZGF0YSA9IGRhdCkgIyBmaXR0aW5nIG3DtCBow6xuaA0Kc3VtbWFyeShnbG0ucG9pcykNCmBgYA0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmxpYnJhcnkoQ09VTlQpDQpQX19kaXNwKGdsbS5wb2lzKQ0KYGBgDQoNCkjDoG0gUF9fZGlzcCBjdW5nIGPhuqVwIGNobyBjaMO6bmcgdGEgUGVhcnNvbiBDaGkyIHN0YXRpc3RpYyB2w6AgUGVhcnNvbiBkaXNwZXJzaW9uIHN0YXRpc3RpYy4gQ8OzIHRo4buDIHTDrW5oIHRo4bunIGPDtG5nIFBlYXJzb24gZGlzcGVyc2lvbiBzdGF0aXN0aWMgYuG6sW5nIGPDoWNoIGzhuqV5IHThu5VuZyBiw6xuaCBwaMawxqFuZyBj4bunYSBwaOG6p24gZMawIFBlYXJzb24gY2hpYSBjaG8gYuG6rWMgdOG7sSBkbyBj4bunYSBtw7QgaMOsbmggKE4gcXVhbiBzw6F0IC0gc+G7kSBiaeG6v24gdHJvbmcgbcO0IGjDrG5oIGJhbyBn4buTbSBpbnRlcmNlcHQpLg0KDQpgYGB7cn0NCnBvaS5kZiA8LSBkZi5yZXNpZHVhbChnbG0ucG9pcykgIyAyMzINClBlYXIuZGlzcC5zdGEgPC0gc3VtKHJlc2lkKGdsbS5wb2lzLCB0eXBlID0gInBlYXJzb24iKV4yKS9wb2kuZGYNClBlYXIuZGlzcC5zdGENCmBgYA0KDQpHacOhIHRy4buLIFBlYXJzb24gZGlzcGVyc2lvbiBzdGF0aXN0aWMgPiAxLCBjw7MgbmdoxKlhIGzDoCBk4buvIGxp4buHdSBn4bq3cCB24bqlbiDEkeG7gSBvdmVyZGlzcGVyc2lvbi4gWGluIGPDsyBt4buZdCBsxrB1IMO9IHbhu5tpIGPDoWMgYuG6oW4gbMOgIGtoaSBk4buvIGxp4buHdSBi4buLIG92ZXJkaXNwZXJzaW9uLCBt4buZdCBnaeG6o2kgcGjDoXAgdGhheSB0aOG6vyBjaG8gbcO0IGjDrG5oIHBvaXNzb24gbMOgIG3DtCBow6xuaCBxdWFzaS1wb2lzc29uLg0KDQpgYGB7cn0NCmdsbS5xdWFzaSA8LSAgZ2xtKHNlaSB+IGxiYXNlICsgdHJ0ICsgdmlzaXQgLCBmYW1pbHkgPSAicXVhc2lwb2lzc29uIiwgZGF0YSA9IGRhdCkNCnN1bW1hcnkoZ2xtLnF1YXNpKSAlPiUgY29lZigpDQpgYGANCg0KSOG7hyBz4buRIGjhu5NpIHF1eSBj4bunYSBtw7QgaMOsbmggZ2xtLnBvaSB2w6AgZ2xtLnF1YXNpIGzDoCBnaeG7kW5nIG5oYXUsIHN0YW5kYXJkIGVycm9yIChTRSkgY+G7p2EgY8OhYyBo4buHIHPhu5EgaOG7k2kgcXV5IOG7nyBtw7QgaMOsbmggZ2xtLnBvaXMgbOG7m24gaMahbiBzbyB24bubaSBTRSDhu58gbcO0IGjDrG5oIGdsbS5wb2lzLg0KDQpDaMO6bmcgdGEgY8OzIHRo4buDIHTDrW5oIHRvw6FuIHPhu7Egc2FpIGtow6FjIFNFIGdp4buvYSAyIG3DtCBow6xuaCBk4buxYSB2w6BvIFBlYXJzb24gZGlzcGVyc2lvbiBzdGF0aXN0aWMuIEtoaSBuaMOibiBTRSBj4bunYSBjw6FjIGjhu4cgc+G7kSBo4buTaSBxdXkg4bufIG3DtCBow6xuaCBwb2lzc29uIHbhu5tpIGPEg24gYuG6rWMgaGFpIGPhu6dhIFBlYXJzb24gZGlzcGVyc2lvbiBzdGF0aXN0aWMgY2jDum5nIHRhIHPhur0gY8OzIMSRxrDhu6NjIFNFIHTGsMahbmcg4bupbmcgdHJvbmcgbcO0IGjDrG5oIHF1YXNpLXBvaXNzb24uDQoNCk5o4buvbmcgZ2nhuqNpIHBow6FwIHRoYXkgdGjhur8ga2jDoWMgY8OzIHRo4buDIGzDoCA6IHBow6JuIHBo4buRaSBQb2lzc29uIEludmVyc2VkIEdhdXNzaWFuIChQSUcpLCBob+G6t2MgTmVnYXRpdmUgYmlub21pYWwgdHlwZSBJSSAoTmjhu4sgdGjhu6ljIMOibSwgTkJJSSkuIA0KDQpDaMO6bmcgdGEgY8OzIHRo4buDIGtp4buDbSBjaOG7qW5nIHbhu4EgdMOtbmggcGjDuSBo4bujcCBj4bunYSBwaMOibiBwaOG7kWkgbmVnYXRpdmUgYmlub21pYWwgSUkgYuG6sW5nIGPDoWNoIDogxJDhu5FpIGNoaeG6v3UgaGlzdG9ncmFtIGNobyBk4buvIGxp4buHdSB0aOG7sWMgdOG6vyB2w6AgaGlzdG9ncmFtIGNobyBt4buZdCBk4buvIGxp4buHdSBtw7QgcGjhu49uZyB0aGVvIHF1eSBsdeG6rXQgTkJJMiwgduG7m2kgY8O5bmcgZ2nDoSB0cuG7iyB0aGFtIHPhu5EgbXUgdsOgIHNpZ21hLiANCkvhur90IHF14bqjIGNobyB0aOG6pXkgMiDEkeG7kyB0aOG7iyBn4bqnbiBuaMawIHTGsMahbmcgxJHGsMahbmcgOg0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmxpYnJhcnkoZ2FtbHNzLmRpc3QpDQoNCnJOQklJKDIzNkwsIG11ID0gbXUsIHNpZ21hID0gc2lnbWEpJT4lYXNfdGliYmxlKCklPiUNCiAgZ2dwbG90KCkrDQogIGdlb21fZGVuc2l0eV9yaWRnZXMoYWVzKHg9dmFsdWUseT0xLGdyb3VwPTEpLHN0YXQgPSAiYmlubGluZSIsYmlud2lkdGg9MSwNCiAgICAgICAgICAgICAgIGRyYXdfYmFzZWxpbmUgPSBGQUxTRSxhbHBoYT0uOCxmaWxsPSJnb2xkIikrDQogIGdlb21fZGVuc2l0eV9yaWRnZXMoYWVzKHg9ZGF0JHNlaSx5PTEsZ3JvdXA9MSksc3RhdCA9ICJiaW5saW5lIixiaW53aWR0aD0xLA0KICAgICAgICAgICAgICAgICAgICAgIGRyYXdfYmFzZWxpbmUgPSBGQUxTRSxhbHBoYT0uNSxmaWxsPSIjYzUwN2ZmIixjb2w9IiNjNTA3ZmYiKSsNCiAgdGhlbWVfYncoKStzY2FsZV95X2NvbnRpbnVvdXMoIkRlbnNpdHkiKSsNCiAgZ2d0aXRsZSgiU2ltdWxhdGVkICh5ZWxsb3cpIHZzIFJlYWwgZGF0YSAodmlvbGV0KSIpDQpgYGANCg0KTmjhuq1uIHjDqXQgOiBWaeG7h2MgbOG7sWEgY2jhu41uIGjhu40gcGjDom4gcGjhu5FpIGPDsyB0aOG7gyBk4buxYSB2w6BvIG3hu5l0IGhv4bq3YyBj4bqjIDIgbOG6rXAgbHXhuq1uIDogMSkgQuG6o24gY2jhuqV0L8OdIG5naMSpYSBsw6JtIHPDoG5nIHbDoCAyKSBUb8OhbiB0aOG7kW5nIGvDqi4gDQoNClbhu5tpICgxKSwgYuG6oW4gxJFhbmcgbmjDrG4gxJHhuqFpIGzGsOG7o25nLCB24bubaSDEkeG6p3kgxJHhu6cgw70gbmdoxKlhLCDEkcahbiB24buLLCB0aGFuZyDEkW/igKYgdGhlbyBxdWFuIMSRaeG7g20geSBo4buNYywgduG7m2kgKDIpLCBvdXRjb21lIGtow7RuZyBjw7JuIGzDoCDEkeG6oWkgbMaw4bujbmcgbuG7r2EgbcOgIGNo4buJIGzDoCBt4buZdCBiaeG6v24gc+G7kSBuZ+G6q3Ugbmhpw6puLCBnacOhIHRy4buLIGPhu6dhIG7DsyBsw6AgMSBjb24gc+G7kSB2w7QgbmdoxKlhIGNo4buLdSBjaGkgcGjhu5FpIGLhu59pIGPDoWMgaMOgbSBt4bqtdCDEkeG7mSB4w6FjIHN14bqldC4gDQoNClRyw6puIHRo4buxYyB04bq/IGhp4bq/bSBraGkgdGEgY8OzIHRo4buDIGxpw6puIGjhu4cgbeG7mXQgY8OhY2ggdMaw4budbmcgbWluaCBnaeG7r2Egw70gbmdoxKlhIGzDom0gc8OgbmcgdsOgIMSR4buLbmggbmdoxKlhIHjDoWMgc3XhuqV0LCB0aMOtIGThu6UgdOG7iSBs4buHIMSRaeG7gXUgdHLhu4sgdGjDoG5oIGPDtG5nIHRyw6puIHThu5VuZyBz4buRIGLhu4duaCBuaMOibiB24burYSBtYW5nIMO9IG5naMSpYSBsw6JtIHPDoG5nLCB24burYSDEkcO6bmcgdGhlbyDEkeG7i25oIG5naMSpYSBj4bunYSBt4buZdCBwaMOibiBwaOG7kWkgQmlub21pYWwgOyBuaMawbmcgdHJvbmcgdHLGsOG7nW5nIGjhu6NwIG7DoHksIHLhuqV0IGtow7MgxJHhu4MgaMOsbmggZHVuZyBz4buRIGPGoW4gxJHhu5luZyBraW5oIGzDoCBz4buRIGzhuqduIHRo4butIG5naGnhu4dtIEJlcm5vdWlsbGkgdGjhuqV0IGLhuqFpIGNobyDEkeG6v24ga2hpIMSR4bqhdCBs4bqnbiB0aMOgbmggY8O0bmcgdGjhu6kgMS9hbHBoYSAtIHRoZW8gxJHhu4tuaCBuZ2jEqWEgcGjDom4gcGjhu5FpIG5o4buLIHRo4bupYyDDom0gdHlwZSAyIDsgdMawxqFuZyB04buxLCBwaMOibiBwaOG7kWkgUG9pc3NvbiBjw7MgdGjhu4MgxJHGsOG7o2MgbOG7sWEgY2jhu41uIGNo4buJIHbDrCB0YSBraMO0bmcgY8OzIHRow7RuZyB0aW4gbsOgbyBraMOhYywgY2jhu6kgdHLDqm4gdGjhu7FjIHThur8ga2jDtG5nIHBo4bqjaSBi4bqldCBj4bupIMSR4bqhaSBsxrDhu6NuZyBuw6BvIGPFqW5nIHRo4buPYSB04bqldCBj4bqjIGdp4bqjIMSR4buLbmggduG7gSB0w61uaCBuZ+G6q3Ugbmhpw6puLCDEkeG7mWMgbOG6rXAgdsOgIGLhuqV0IGJp4bq/biBj4bunYSBwaMOibiBwaOG7kWkgUG9pc3Nvbi4gVGjDrSBk4bulIDogbcO5YSB4dcOibiBsw6AgdGjhu51pIMSRaeG7g20gYuG7h25oIG5ow6JuIGhlbiBk4buLIOG7qW5nIHbhu5tpIHBo4bqlbiBob2Ega2jhu59pIHBow6F0IG5oaeG7gXUgY8ahbiBoZW4gaMahbiBzbyB24bubaSBtw7lhIGjDqCwgeMOhYyBzdeG6pXQgY+G7p2EgY8ahbiBoZW4gxJHDoyBi4buLIHRoYXkgxJHhu5VpIGNo4bupIGtow7RuZyBjw7JuIGjhurFuZyDEkeG7i25oIG7hu69hLiANCg0KTmjGsCB24bqteSwgcGjDom4gcGjhu5FpIG5o4buLIHRo4bupYyDDom0gxJHGsOG7o2MgY2jhu41uIGNobyBuZ2hpw6puIGPhu6l1IG7DoHkgbmjGsCBt4buZdCBwaMawxqFuZyB0aeG7h24sIMSR4buDIMSRaeG7gXUgY2jhu4luaCB0w61uaCBxdcOhIHBow6JuIHTDoW4gY+G7p2Egc+G7kSBsaeG7h3UsIGNo4bupIGtow7RuZyBwaOG6o2kgbmjGsCBt4buZdCBxdXkgbHXhuq10Lg0KDQojIyBNw7QgaMOsbmggTkJJSSA6IFThu6sgUkVNTCDEkeG6v24gQmF5ZXMNCg0KQ2hvIGLDoGkgdG/DoW4gaOG7k2kgcXV5IHbhu5tpIG3hu6VjIHRpw6p1IMaw4bubYyB0w61uaCBnacOhIHRy4buLIGJp4bq/biBr4bq/dCBxdeG6oyB5ICh0aGF5IGNobyB4KSBjw7MgcGjDom4gcGjhu5FpIE5CSSB24bubaSB0aGFtIHPhu5EgdHJ1bmcgYsOsbmggTXUgdsOgIGtp4buDdSBow6xuaCAxL2FscGhhIChhbHBoYSB0aGF5IGNobyAxL3IsIHRoZW8gxJHhu4tuaCBuZ2jEqWEgduG7gSBz4buRIGzhuqduIHRow6BuaCBjw7RuZyB0aOG7qSByKS4gUGjDom4gcGjhu5FpIE5lZ2F0aXZlIGJpbm9taWFsIGPDsm4gY8OzIHRo4buDIMSRxrDhu6NjIHRyw6xuaCBiw6B5IGTGsOG7m2kgaMOsbmggdGjhu6ljIA0KDQokJFAoeXxcbXUsXGFscGhhKSA9IFxiaW5vbXt5KzEvXGFscGhhIC0xfXsxL1xhbHBoYS0xfSggXGZyYWN7MX17MStcYWxwaGEgXG11fSleezEvXGFscGhhfSAoXGZyYWN7XGFscGhhIFxtdX17MStcYWxwaGEgXG11fSlee3l9JCQNCg0KSGF5IGjDrG5oIHRo4bupYyBo4buXbiBo4bujcCBow6BtIFBvaXNzb24gdsOgIEdhbW1hIDoNCg0KJCRQKHl8XG11LFxhbHBoYSkgPSBcZnJhYyB7XEdhbW1hKHkrMS9cYWxwaGEpfXtcR2FtbWEoeSsxKVxHYW1tYSgxL1xhbHBoYSl9KCBcZnJhY3sxfXsxK1xhbHBoYSBcbXV9KV57MS9cYWxwaGF9IChcZnJhY3tcYWxwaGEgXG11fXsxK1xhbHBoYSBcbXV9KV57eX0kJA0KxJBp4buBdSBraeG7h24gZHV5IG5o4bqldCDhu58gxJHDonkgbMOgIE11IHbDoCBBbHBoYSA+IDAuIE5oxrAgxJHDoyBuw7NpIOG7nyB0csOqbiwgMS9hbHBoYSA9IHIga2jDtG5nIGLhuq90IGJ14buZYyBwaOG6o2kgbMOgIHPhu5Egbmd1ecOqbiAuDQoNCk3DtCBow6xuaCBuZWdhdGl2ZSBiaW5vbWlhbCB0eXBlIDIgY8OybiDEkcaw4bujYyBn4buNaSBtw7QgaMOsbmggcG9pc3Nvbi1nYW1tYSBkaXN0cmlidXRpb24gbWl4dHVyZSBtb2RlbC4gVuG7m2kgYmnhur9uIHkgdHXDom4gdGhlbyBwaMOibiBwaOG7kWkgTkIyIHRow6wgY2jDum5nIHRhIGPhuqduIGzGsHUgw70gbmjhu69uZyDEkWnhu4NtIHNhdToNCg0KKyBHaeG6oyDEkeG7i25oIHkgdHXDom4gdGhlbyBwaMOibiBwaOG7kWkgcG9pc3NvbiB24bubaSBtZWFuID0gdmFyaWFuY2UgPSBNdQ0KDQorIE11IHPhur0gdHXDom4gdGhlbyBwaMOibiBwaOG7kWkgZ2FtbWEgduG7m2kgTWVhbiA9IE11IHbDoA0KDQokJHZhcmlhbmNlID0gXG11ICsgXGZyYWN7XG11XjJ9e3Z9JCQNCg0K4bueIMSRw6J5IHYgxJHGsOG7o2MgZ+G7jWkgbMOgIHNoYXBlIHBhcmFtZXRlciAodGhhbSBz4buRIGtp4buDdSBow6xuaCkgdHJvbmcgcGjDom4gcGjhu5FpIGdhbW1hLCBuw7MgdMawxqFuZyDhu6luZyB24bubaSB0aGFtIHPhu5EgciB0cm9uZyBwaMOibiBwaOG7kWkgTkJJIGNow61uaCB04bqvYyAoY8OybiBn4buNaSBsw6AgaW5kaXJlY3QgZGlzcGVyc2lvbiBwYXJhbWV0ZXIgdHJvbmcgcGjDom4gcGjhu5FpIE5CKTsgdGEgdGhheSByIGhheSB2IGLhurFuZyBhbHBoYSBsw6AgZGlyZWN0IGRpc3BlcnNpb24gcGFyYW1ldGVyLiANCg0KTmjGsCB24bqteSDEkeG7kWkgduG7m2kgbcO0IGjDrG5oIE5CLCBjw7MgMiB0aGFtIHPhu5EgY+G6p24gxJHGsOG7o2MgxrDhu5tjIHTDrW5oIGzDoCBNdSAodHJ1bmcgYsOsbmgpIHbDoCBkaXNwZXJzaW9uIHBhcmFtZXRlciAodiBob+G6t2MgYWxwaGEpLiBMxrB1IMO9IHLhurFuZyBraGkgc+G7rSBk4bulbmcgUiwgdGjDrCB2IMSRxrDhu6NjIGfhu41pIGLhurFuZyB0w6puIHRoYXkgdGjhur8gbMOgIFRoZXRhIOG7nyB0cm9uZyBwaOG6p24ga+G6v3QgcXXhuqMgY+G7p2EgbcO0IGjDrG5oLiBN4buZdCBz4buRIHBo4bqnbiBt4buBbSBraMOhYyBuaMawIFN0YXRhLCBTQVMga2jDtG5nIHTDrW5oIHRvw6FuIHRoZXRhIGhheSB2IG3DoCBz4bq9IHTDrW5oIGFscGhhIHbDoCBn4buNaSDEkcOzIGzDoCBkaXNwZXJzaW9uIHBhcmFtZXRlciB0cm9uZyBwaOG6p24ga+G6v3QgcXXhuqMuIA0KDQpLaGkgTkIyIG1vZGVsIMSRxrDhu6NjIG3DtCBow6xuaCBow7NhIGTGsOG7m2kgZOG6oW5nIGvhur90IGjhu6NwIGPhu6dhIGdhbW1hLXBvaXNzb24gdGjDrCB2YXJpYW5jZSBj4bunYSBtw7QgaMOsbmggTkIgc+G6vSBi4bqxbmcgdOG7lW5nIHZhcmlhbmNlIGPhu6dhIHBow6JuIHBo4buRaSBwb2lzc29uIHbDoCBwaMOibiBwaOG7kWkgZ2FtbWEuDQoNClRoYW0gc+G7kSBNdSBjw7MgdGjhu4MgxJHGsOG7o2MgxrDhu5tjIGzGsOG7o25nIGLhurFuZyBtw7QgaMOsbmggR0xNIHRow7RuZyBxdWEgMSBsaW5rIGZ1bmN0aW9uIGzDoCBow6BtIGxvZ2FyaXQgdOG7sSBuaGnDqm4gKExOKSA6DQoNCiQkbG4oXG11KSBcc2ltIFxiZXRhX3swfSArIFxzdW1fe2k9MX1ee0t9XGJldGFfe2l9IHhfe2l9JCQNCg0KVuG7m2kgYmV0YTAgbMOgIGludGVyY2VwdCB2w6AgYmV0YShpKSBsw6AgdGhhbSBz4buRIGjhu5NpIHF1eSBjaG8gYmnhur9uIHPhu5EgxJHhu5ljIGzhuq1wIHgoaSkgdHJvbmcgbW9kZWwgbWF0cml4IFggZ+G7k20gayBiaeG6v24gc+G7kQ0KDQpUaGVvIHRyxrDhu51uZyBwaMOhaSDGsOG7m2MgdMOtbmggaOG7o3AgbMO9IGPhu7FjIMSR4bqhaSAoTWF4aW11bSBsaWtlbGlob29kKSwgdGEgY8OzIHRo4buDIMaw4bubYyBsxrDhu6NuZyB0aGFtIHPhu5EgYWxwaGEgdsOgIGJldGEgKGNobyBtdSkgdGjDtG5nIHF1YSBow6BtIGxpa2VsaWhvb2QgOg0KDQokJEwoXGJldGF8XGFscGhhICxcYmV0YSApPSBccHJvZF97aT0xfV57bn1QKHlfe2l9KSA9IFxwcm9kX3tpPTF9XntufSBcZnJhYyB7XEdhbW1hKHlfe2l9KzEvXGFscGhhKX17XEdhbW1hKHlfe2l9KzEpXEdhbW1hKDEvXGFscGhhKX0oIFxmcmFjezF9ezErXGFscGhhIGVee3hfe2l9XGJldGF9fSleezEvXGFscGhhfSAoXGZyYWN7XGFscGhhIGVee3hfe2l9XGJldGF9fXsxK1xhbHBoYSBlXnt4X3tpfVxiZXRhfX0pXnt5X3tpfX0kJA0KDQpHacOhIHRy4buLIGFscGhhIHbDoCBiZXRhIGjhu6NwIGzDvSBuaOG6pXQgdMawxqFuZyDhu6luZyB24bubaSBMIMSR4bqhdCBj4buxYyDEkeG6oWkuDQoNCk3DtCBow6xuaCBtaXhlZCBOQiB24bubaSByYW5kb20gZWZmZWN0IHRoZW8gdHLGsOG7nW5nIHBow6FpIGZyZXF1ZW50aXN0IGPDsyB0aOG7gyDEkcaw4bujYyBmaXQgduG7m2kgcGFja2FnZSBnYW1sc3MgaGF5IGxtZTQuIMSQ4buDIMSRxqFuIGdp4bqjbiB24buBIG3hurd0IGPDuiBwaMOhcCwgY2jDum5nIHRhIHPhur0gZMO5bmcgcGFja2FnZSBsbWU0LiBW4bubaSBi4buZIGThu68gbGnhu4d1IGVwaWwsIGNow7puZyB0YSBz4bq9IHhlbSB4w6l0IOG6o25oIGjGsOG7n25nIGPhu6dhIHRydCwgbGJhc2UgdsOgIHZpc2l0IGzDqm4gc+G7kSBs4bqnbiBi4buLIMSR4buZbmcga2luaCB24bubaSByYW5kb20gaW50ZXJjZXB0IGNobyB04burbmcgYuG7h25oIG5ow6JuLg0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCnJlcXVpcmUobG1lNCkNCg0KbmIuYWRtYiA8LSBnbG1lci5uYihzZWkgfiB0cnQrbGJhc2UrdmlzaXQgKyAoMXxzdWJqZWN0KSwgZGF0YSA9IGRhdCApDQpzdW1tYXJ5KG5iLmFkbWIpDQpgYGANCg0KS2hpIGThu7FuZyBtw7QgaMOsbmggaOG7k2kgcXV5IG5o4buLIHRo4bupYyDDom0gdGhlbyBCYXllcywgbeG7pWMgdGnDqnUgY+G7p2EgY2jDum5nIHRhIGzDoCDGsOG7m2MgdMOtbmggcGjDom4gcGjhu5FpIGjhuq11IMSR4buLbmggY2hvIHRoYW0gc+G7kSBiZXRhIChjaG8gcGjDqXAgxrDhu5tjIHTDrW5oIE11KSB2w6AgdGhhbSBz4buRIGFscGhhIGzDoCBraeG7g3UgaMOsbmguIE3DtCBow6xuaCBuw6B5IGPFqW5nIGPDsyB0aOG7gyBjaOG7qWEgcmFuZG9tIGVmZmVjdC4gUHJpb3IgdsOgIGxpa2VsaWhvb2QgY+G7p2EgbcO0IGjDrG5oIEJheWVzIE5CSSDEkcaw4bujYyBtaW5oIGjhu41hIHRyb25nIGjDrG5oIHNhdToNCg0KIVtdKE5CSUJheWVzMi5wbmcpDQoNCk3DtCBow6xuaCBtaXhlZCBOQiB24bubaSByYW5kb20gaW50ZXJjZXB0IGNow7puZyB0YSBj4bqnbiDGsOG7m2MgdMOtbmggY8OzIGThuqFuZzoNCg0KJCRcbXVfe2lqfSA9IGVeXGxlZnQgKCBcYmV0YV97MH0gKyBcYmV0YV97dHJ0fSAgeF97dHJ0fSArIFxiZXRhX3tsYmFzZX0gIHhfe2xiYXNlfSArIFxiZXRhX3t2aXNpdH0gIHhfe3Zpc2l0fSArIHN1Yl9pICsgXHZhcmVwc2lsb25fe2lqfSAgXHJpZ2h0ICkkJA0KDQojIEzhu7FhIGNo4buNbiBtw7QgaMOsbmggTkJJMiB0cm9uZyBTVEFODQoNCkLDonkgZ2nhu50gY2jDum5nIHRhIHPhur0gc+G7rSBk4bulbmcgU3RhbiDEkeG7gyDEkeG7gyBmaXQgbcO0IGjDrG5oIHRow7RuZyBxdWEgcGFja2FnZSByc3Rhbi4NCg0KKyBCxrDhu5tjIDE6IFZp4bq/dCBTdGFuIGNvZGUuDQoNCuG7niBuaOG7r25nIGLDoGkgdGjhu7FjIGjDoG5oIHRyxrDhu5tjIGPDoWMgYuG6oW4gxJHDoyBsw6BtIHF1ZW4gduG7m2kgY8OhY2ggdmnhur90IHN0YW4gY29kZSB0aGVvIHThu6tuZyBibG9jaywgbsOqbiBjaMO6bmcgdMO0aSB4aW4gcGjDqXAga2jDtG5nIG5o4bqvYyBs4bqhaSBs4bqnbiBu4buvYS4g4bueIMSRw6J5IGNow7puZyB0w7RpIHPhur0gZ2hpIGNow7ogdHLhu7FjIHRp4bq/cCB0cm9uZyB04burbmcgYmxvY2sgxJHhu4MgZ2nhuqNpIHRow61jaCDDvSBuZ2jEqWEgY+G7p2EgdOG7q25nIGTDsm5nIGzhu4duaC4NCg0KYGBge3IsbWVzc2FnZSA9IEZBTFNFLHdhcm5pbmc9RkFMU0V9DQpzdGFuX05CSS5zdHJpID0iDQpkYXRhIHsNCmludDxsb3dlcj0wPiBOOyAgICAgICAgLy9ubyBvZiBvYnMNCmludDxsb3dlcj0wPiBLOyAgICAgICAgLy9ubyBvZiBwYXJhbWV0ZXJzIChpbnRlcmNlcHQrcHJlZGljdG9ycykNCmludDxsb3dlcj0wPiBOc3ViOyAgICAgLy9ubyBvZiBzdWJzDQppbnQ8bG93ZXI9MD4gc3ViW05dOyAgIC8vIGluZGV4IGZvciBzdWJqZWN0DQppbnQ8bG93ZXI9MD4geVtOXTsgICAgIC8vZGVwZW5kZW50IHZhcmlhYmxlDQptYXRyaXhbTixLXSBYOyAgICAgICAgIC8vbW9kZWwgbWF0cml4DQp2ZWN0b3JbS10gYnByaTsgICAgICAgIA0KbWF0cml4W0ssS10gQnByaTsNCnZlY3RvcltOc3ViXSBhcHJpOw0KbWF0cml4W05zdWIsIE5zdWJdIEFwcmk7DQp9DQpwYXJhbWV0ZXJzIHsNCnZlY3RvcltLXSBiZXRhOyAgICAgICAgICAgIC8vYmV0YVsxXSBpcyBwb3B1bGF0aW9uIGludGVyY2VwdA0KdmVjdG9yW05zdWJdIHJlX2l0Y2VwdDsgICAgLy8gcmFuZG9tIGludGVyY2VwdA0KcmVhbDxsb3dlcj0wPiBzaWdtYV9yZTsgICAgICAgICAgICAgLy9TRCByYW5kb20gaW50ZXJjZXB0DQpyZWFsPGxvd2VyPTA+IGFscGhhOyAgICAgICAgICAgICAgICAvL2RpcmVjdCBkaXNwZXJzaW9uIHBhcmFtZXRlcg0KfQ0KdHJhbnNmb3JtZWQgcGFyYW1ldGVycyB7DQp2ZWN0b3JbTl0gbG9nX2xpazsNCnZlY3RvcltOXSBldGE7ICAgICAgICAgICAgIC8vbGluZWFyIHByZWRpY3Rvcg0KdmVjdG9yW05dIG11OyAgICAgICAgICAgICAgLy9lc3RpbWF0ZWQgbWVhbiBvZiBkZXBlbmRlbnQgdmFyaWFibGUNCmV0YSA9IFgqYmV0YTsNCmZvcihpIGluIDE6Tikgew0KbXVbaV0gPSBleHAoZXRhW2ldICsgcmVfaXRjZXB0W3N1YltpXV0pOw0KbG9nX2xpa1tpXSA9IG5lZ19iaW5vbWlhbF8yX2xwbWYoeVtpXXxtdVtpXSwxL2FscGhhKTsNCn0NCn0NCm1vZGVsIHsNCi8vcHJpb3INCmFscGhhIH4gY2F1Y2h5KDAsMjUpOw0Kc2lnbWFfcmUgfiBjYXVjaHkoMCw1KTsNCmJldGEgfiBtdWx0aV9ub3JtYWwoYnByaSwgQnByaSk7ICAgICAvL21lYW4gPSAwOyBzZCA9IDENCnJlX2l0Y2VwdCB+IG11bHRpX25vcm1hbChhcHJpLCBBcHJpKnNpZ21hX3JlKTsgICAvL21lYW4gPSAwDQovL2xpa2VsaWhvb2QNCnkgfiBuZWdfYmlub21pYWxfMihtdSwgMS9hbHBoYSk7IC8vdGhldGEgPSAxL2FscGhhDQp9DQoiDQpgYGANCg0KKyBCxrDhu5tjIDI6IEzhu7FhIGNo4buNbiBtw7QgaMOsbmggYuG6sW5nIHbDsm5nIGzhurdwIFdBSUMNCg0KVHLGsOG7m2Mga2hpIHRo4buxYyBz4buxIGThu7FuZyBtw7QgaMOsbmgsIGNow7puZyB0w7RpIHPhur0gZ2nhu5tpIHRoaeG7h3UgduG7gSBwaMawxqFuZyBwaMOhcCBs4buxYSBjaOG7jW4gbcO0IGjDrG5oIChiaeG6v24gc+G7kSkgZOG7sWEgdsOgbyB0acOqdSBjaMOtIFdBSUMgKFdhdGFuYWJlLUlrYWlrZSBpbmZvcm1hdGlvbiBjcml0ZXJpYSkgdsOgIExPT0lDIChMZWF2ZS1vbmUtb3V0IGluZm9ybWF0aW9uIGNyaXRlcmlhKS4gQ8ahIGNo4bq/IMSR4bqxbmcgc2F1IG5o4buvbmcgdHLhu4sgc+G7kSBuw6B5IHPhur0gxJHGsOG7o2MgZ2nhuqNpIHRow61jaCByw7UgaMahbiB0cm9uZyBuaOG7r25nIGLDoGkgdGnhur9wIHRoZW8uDQoNClThuqFpIHRo4budaSDEkWnhu4NtIG7DoHksIGNow7puZyB0YSB04bqhbSBjaOG6pXAgbmjhuq1uIHF1eSB04bqvYyA6IFdBSUMgdsOgIExPT0lDIMSRxrDhu6NjIHTDrW5oIHThu6sgTG9nLWxpa2VsaWhvb2QgY+G7p2EgbcO0IGjDrG5oIEJheWVzLCB2w6AgbcO0IGjDrG5oIHBow7kgaOG7o3AgbmjhuqV0IHbhu5tpIGThu68gbGnhu4d1IHPhur0gY8OzIFdBSUMgdsOgIExPT0lDIHRo4bqlcCBuaOG6pXQuIGLhurFuZyBjw6FjaCBzbyBzw6FuaCBXQUlDIHbDoCBMT09JQyBnaeG7r2EgY8OhYyBwaGnDqm4gYuG6o24gbW9kZWwga2jDoWMgbmhhdSwgdGEgY8OzIHRo4buDIGNo4buNbiDEkcaw4bujYyBtb2RlbCBo4bujcCBsw70gbmjhuqV0LiBDw6FjaCBsw6BtIG7DoHkgZ2nhu5FuZyBuaMawIHF1eSB0csOsbmggU3RlcC13aXNlIGThu7FhIHbDoG8gQUlDICBoYXkgQklDIHRyb25nIHRyxrDhu51uZyBwaMOhaSBNYXhpbXVtIExpa2VsaWhvb2QgbcOgIGLhuqFuIMSRw6MgYmnhur90Lg0KDQpDaMO6bmcgdMO0aSB2aeG6v3QgMSB2w7JuZyBs4bq3cCwgbOG6p24gbMaw4bujdCBk4buxbmcgNyBtw7QgaMOsbmggTkJJIEJheWVzIGNobyB04bqldCBj4bqjIHThu5UgaOG7o3AgY8OzIHRo4buDIGdp4buvYSAzIGJp4bq/biBz4buROiBUcmVhdG1lbnQsIExvZ2Jhc2UgdsOgIFZpc2l0LCB2w6AgdMOtbmggV0FJQywgTE9PSUMgY2hvIG3hu5dpIG3DtCBow6xuaC4NCg0KYGBge3J9DQpsaWJyYXJ5KGxvbykNCmBgYA0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmxpYnJhcnkocnN0YW4pDQoNCmxpYnJhcnkocGFyYWxsZWwpDQoNCmNvbSA8LSBjKGNvbWJuKDI6NCwgMSwgc2ltcGxpZnk9RkFMU0UpLCANCiAgICAgICAgIGNvbWJuKDI6NCwgMiwgc2ltcGxpZnk9RkFMU0UpLCANCiAgICAgICAgIGNvbWJuKDI6NCwgMywgc2ltcGxpZnk9RkFMU0UpKQ0KDQp0cnQ9ZGF0JHRydA0KbGJhc2U9ZGF0JGxiYXNlDQp2aXNpdD1kYXQkdmlzaXQNCg0KWG1hdCA9IG1vZGVsLm1hdHJpeCh+IHRydCtsYmFzZSt2aXNpdCkNCg0Kb3V0cHV0PWRhdGEuZnJhbWUoDQogIFdBSUM9cmVwKE5BLDcpLA0KICBTRVdBSUM9cmVwKE5BLDcpLA0KICBMT09JQz1yZXAoTkEsNyksDQogIFNFTE9PSUM9cmVwKE5BLDcpLA0KICBNb2RlbElEPWMoMTo3KQ0KKQ0KDQpOc3ViID0gbmxldmVscyhhcy5mYWN0b3IoZGF0JHN1YmplY3QpKSAgIyBubyBvZiBzdWJqZWN0cw0KTj1ucm93KGRhdCkNCg0KUGJhciA8LSB3aW5Qcm9ncmVzc0Jhcih0aXRsZSA9ICJXQUlDIGxvb3AiLCBtaW4gPSAwLA0KICAgICAgICAgICAgICAgICAgICAgICBtYXggPSBsZW5ndGgoY29tKSwgd2lkdGggPSAzNTApDQoNCmZvcihtIGluIDE6bGVuZ3RoKGNvbSkpIHsNCiAgWD0gWG1hdFssYygxLCB1bmxpc3QoY29tW1ttXV0pKV0NCiAgSz0gbmNvbChYKQ0KICB0ZW1wZGF0IDwtIGxpc3QoeSA9IGRhdCRzZWksICAjICByZXNwb25zZSB2YXJpYWJsZQ0KICAgICAgICAgICAgICAgICAgICAgICAgTj1OLCAgICAgICAgIyBubyBvZiBvYnMNCiAgICAgICAgICAgICAgICAgICAgICAgIFg9WCwgIyBtb2RlbCBtYXRyaXggb2YgcHJlZGljdG9ycyBhbmQgaW50ZXJjZXB0DQogICAgICAgICAgICAgICAgICAgICAgICBLPUssICAgICAgICAgIyBubyBvZiBwcmVkaWN0b3JzICsgaW50ZXJjZXB0DQogICAgICAgICAgICAgICAgICAgICAgICBOc3ViID0gTnN1YiwgICAgICAgICAgICAjIG5vIG9mIHN1YnMgKCA1OSA9IG5vIG9mIHJhbmRvbSBlZmZlY3QgVkFMVUVzKQ0KICAgICAgICAgICAgICAgICAgICAgICAgc3ViID0gZGF0JHN1YmplY3QlPiVhcy5udW1lcmljKCksICAgICAgIyBTdWJqZWN0IElEDQogICAgICAgICAgICAgICAgICAgICAgICBicHJpID0gcmVwKDAsIEspLCAgICAgICAjIFRvIGRlZmluZSBwcmlvciBmb3IgYmV0YQ0KICAgICAgICAgICAgICAgICAgICAgICAgQnByaSA9IGRpYWcoMSwgSyksICAgICAgIyBQcmlvciBtYXRyaXggb2YgYmV0YSANCiAgICAgICAgICAgICAgICAgICAgICAgIGFwcmkgPSByZXAoMCwgTnN1YiksICAgICMgVG8gZGVmaW5lIHByaW9yIGZvciBpbnRlcmNlcHQNCiAgICAgICAgICAgICAgICAgICAgICAgIEFwcmkgPSBkaWFnKDEsIE5zdWIgKSkgIyBQcmlvciBtYXRyaXggb2YgaW50ZXJjZXB0DQoNCnNldC5zZWVkKDEyMykNCg0KdGVtcG1vZCA8LSBzdGFuKG1vZGVsX2NvZGUgPSBzdGFuX05CSS5zdHJpLCBkYXRhID0gdGVtcGRhdCwgDQogICAgICAgICAgICAgICAgaXRlciA9IDE1MDAsIHdhcm11cCA9IDUwMCkgDQoNCncxPXRlbXBtb2QlPiVleHRyYWN0X2xvZ19saWsoLixwYXJhbWV0ZXJfbmFtZSA9ImxvZ19saWsiKSU+JXdhaWMoKSANCmwxPXRlbXBtb2QlPiVleHRyYWN0X2xvZ19saWsoLixwYXJhbWV0ZXJfbmFtZSA9ImxvZ19saWsiKSU+JWxvbygpIA0KDQpvdXRwdXQkV0FJQ1ttXTwtdzEkd2FpYw0Kb3V0cHV0JFNFV0FJQ1ttXTwtdzEkc2Vfd2FpYw0Kb3V0cHV0JExPT0lDW21dPC1sMSRsb29pYw0Kb3V0cHV0JFNFTE9PSUNbbV08LWwxJHNlX2xvb2ljDQoNCnNldFdpblByb2dyZXNzQmFyKFBiYXIsIG0sIHRpdGxlPXBhc3RlKCByb3VuZChtL2xlbmd0aChjb20pKjEwMCwgMCksDQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIiUgSG/DoG4gdGjDoG5oIikpDQp9DQoNCmxpYnJhcnkodmlyaWRpcykNCg0Kb3V0cHV0JT4lZ2dwbG90KGFlcyh4PXJlb3JkZXIoTW9kZWxJRCxXQUlDKSx5PVdBSUMpKSsNCiAgZ2VvbV9wb2ludChzaGFwZT0yMSxhZXMoZmlsbD1mYWN0b3IoTW9kZWxJRCksDQogICAgICAgICAgICAgICAgICAgICAgICAgIHNpemU9V0FJQyxjb2w9ZmFjdG9yKE1vZGVsSUQpLA0KICAgICAgICAgICAgICAgICAgICAgICAgICBzdHJva2U9U0VXQUlDLzEuNSksDQogICAgICAgICAgICAgc2hvdy5sZWdlbmQgPSBGLGFscGhhPTAuNSkrDQogIGdlb21fcG9pbnQoc2hhcGU9MjEsYWVzKGZpbGw9ZmFjdG9yKE1vZGVsSUQpLA0KICAgICAgICAgICAgICAgICAgICAgICAgICBzaXplPVdBSUMpLA0KICAgICAgICAgICAgIHNob3cubGVnZW5kID0gRikrDQogIHRoZW1lX2J3KCkrc2NhbGVfeF9kaXNjcmV0ZSgiTW9kZWwgdmVyc2lvbiIpKw0KICBzY2FsZV9jb2xvdXJfdmlyaWRpcyhvcHRpb249IkMiLGRpc2NyZXRlPVQpKw0KICBzY2FsZV9maWxsX3ZpcmlkaXMob3B0aW9uPSJDIixkaXNjcmV0ZSA9IFQpDQogIA0Kb3V0cHV0JT4lZ2dwbG90KGFlcyh4PXJlb3JkZXIoTW9kZWxJRCxMT09JQykseT1MT09JQykpKw0KICBnZW9tX3BvaW50KHNoYXBlPTIxLGFlcyhmaWxsPWZhY3RvcihNb2RlbElEKSwNCiAgICAgICAgICAgICAgICAgICAgICAgICAgc2l6ZT1MT09JQyxjb2w9ZmFjdG9yKE1vZGVsSUQpLA0KICAgICAgICAgICAgICAgICAgICAgICAgICBzdHJva2U9U0VMT09JQy8xLjUpLA0KICAgICAgICAgICAgIHNob3cubGVnZW5kID0gRixhbHBoYT0wLjUpKw0KICBnZW9tX3BvaW50KHNoYXBlPTIxLGFlcyhmaWxsPWZhY3RvcihNb2RlbElEKSwNCiAgICAgICAgICAgICAgICAgICAgICAgICAgc2l6ZT1MT09JQyksDQogICAgICAgICAgICAgc2hvdy5sZWdlbmQgPSBGKSsNCiAgdGhlbWVfYncoKStzY2FsZV94X2Rpc2NyZXRlKCJNb2RlbCB2ZXJzaW9uIikrDQogIHNjYWxlX2NvbG91cl92aXJpZGlzKG9wdGlvbj0iQyIsZGlzY3JldGU9VCkrDQogIHNjYWxlX2ZpbGxfdmlyaWRpcyhvcHRpb249IkMiLGRpc2NyZXRlID0gVCkNCg0KYGBgDQoNCkvhur90IHF14bqjIGNobyB0aOG6pXk6IE3DtCBow6xuaCBwaMO5IGjhu6NwIG5o4bqldCBsw6AgbcO0IGjDrG5oIHPhu5EgNywgbuG7mWkgZHVuZyBj4bunYSBuw7MgbMOgOg0KDQpgYGB7cn0NCmNvbG5hbWVzKFgpDQpgYGANCg0KTmjGsCB24bqteSB0YSBjw7MgdGjhu4MgY2jhuq9jIGNo4bqvbiB0cm9uZyBtw7QgaMOsbmggbsOgeSBjw7MgMyBiaeG6v24gc+G7kSBuaMawIGdp4bqjIHRodXnhur90IGJhbiDEkeG6p3UNCg0KVHLGsOG7m2Mga2hpIGThu7FuZyBtw7QgaMOsbmgsIGNow7puZyB0YSBjw7MgdGjhu4MgdGjEg20gZMOyIGThu68gbGnhu4d1IDEgY2jDunQ6DQoNCisgVMawxqFuZyBwaOG6o24gZ2nhu69hIDIgcGjDom4gbmjDs20gdHLhu4sgbGnhu4d1IHThuqFpIDQgdGjhu51pIMSRaeG7g20gKHZpc2l0IGtow6FjIG5oYXU6KQ0KDQpgYGB7cn0NCmRlbi5yaWQgPC0gZGF0ICU+JSBnZ3Bsb3QoLiwgYWVzKHggPSBzZWksIHkgPSB0cnQsIGZpbGw9IHRydCkpICsgDQogIGdlb21fZGVuc2l0eV9yaWRnZXMoc3RhdCA9ICJiaW5saW5lIiwgc2NhbGUgPSAxLCBiaW53aWR0aD0yLCBkcmF3X2Jhc2VsaW5lID0gRkFMU0UpICsNCiAgbGFicyh4PSJDb3VudCIsIHkgPSAiIikgKyBzY2FsZV9maWxsX2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IikgKw0KICBjb29yZF9mbGlwKCkgKw0KICBteV90aGVtZSgxMCkgKw0KICBnZW9tX3ZsaW5lKHhpbnRlcmNlcHQgPSBtZWRpYW4oZGF0JGJhc2UpLGxpbmV0eXBlPTIsY29sPSJibHVlIikrDQogIHRoZW1lKGF4aXMudGV4dC54ID0gZWxlbWVudF90ZXh0KGFuZ2xlID0gNDUsIHZqdXN0ID0gMSwgaGp1c3Q9MSkpICsNCiAgdGhlbWUoYXhpcy50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSAxMiwgY29sb3IgPSAiYmxhY2siKSwNCiAgYXhpcy50aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gMTUsIGNvbG9yID0gImJsYWNrIikpICsNCiAgZmFjZXRfd3JhcCh+IHBlcmlvZCwgbmNvbD00KQ0KDQpkZW4ucmlkIA0KYGBgDQoNClRhIGPDsyB0aOG7gyB0aOG6pXkgcuG6sW5nIHTDrG5oIHRy4bqhbmcgYuG7h25oIGzDvSB0csaw4bubYyBraGkgxJFp4buBdSB0cuG7iyAobG9nIGJhc2VsaW5lKWPDsyDhuqNuaCBoxrDhu59uZyDDrXQgbmhp4buBdSDEkeG6v24gaGnhu4d1IHF14bqjIMSRaeG7gXUgdHLhu4sgY+G7p2EgcHJvZ2FiaWRlOg0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmcuc21vb3RoIDwtIGRhdCU+JWdncGxvdChhZXMoeD10cnQseT1zcXJ0KHNlaSksIGdyb3VwPWxiYXNlKSkrDQogIGdlb21fcG9pbnQoYWxwaGE9MC44LGFlcyhjb2xvcj1sYmFzZSkpKw0KICBnZW9tX3Ntb290aChzaXplPTEsYWVzKGNvbG9yPWxiYXNlKSxzZT1GLG1ldGhvZD0ibG0iLHNob3cubGVnZW5kID0gRikrDQogIG15X3RoZW1lKDEyKStzY2FsZV9jb2xvcl92aXJpZGlzKG9wdGlvbj0iQSIsIGRpcmVjdGlvbiA9IDEpKw0KICB0aGVtZShheGlzLnRleHQgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDExLCBjb2xvciA9ICJibGFjayIpLA0KICBheGlzLnRpdGxlID0gZWxlbWVudF90ZXh0KHNpemUgPSAxNSwgY29sb3IgPSAiYmxhY2siKSkgKw0KICBmYWNldF93cmFwKH4gcGVyaW9kLG5jb2w9NCkNCmcuc21vb3RoIA0KYGBgDQoNClRhIGPDsyB0aOG7gyBraOG6o28gc8OhdCBoaeG7h3Ug4bupbmcgY+G7p2EgVGjhu51pIGdpYW4gbMOqbiBMb2dhcml0IGPhu6dhIGJp4bq/biBz4buRIGvhur90IHF14bqjLCBjw7MgduG6uyBuaMawIGhp4buHdSDhu6luZyBuw6B5IHLhuqV0IHnhur91Og0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmcudHJ0IDwtIGRhdCAlPiUgZ2dwbG90KC4sYWVzKHg9cGVyaW9kLCB5PWxvZyhzZWkgKyAwLjAxKSwgZ3JvdXAgPSB0cnQpKSArIA0KICAgICBnZW9tX3BvaW50KGFlcyhjb2xvciA9IHRydCksIHNpemUgPSAyICkgKw0KICAgICBnZW9tX3BvaW50KHNoYXBlID0xLCBzaXplID0gMiwgY29sb3IgPSAiYmxhY2siKSArDQogICAgIGdlb21fc21vb3RoKGFlcyhjb2xvcj0gdHJ0LCBmaWxsID0gdHJ0KSwgc2l6ZSA9IDEsIGFscGhhPSAwLjMsIHNlID0gVCwgbWV0aG9kPSJsbSIpICsNCiAgICAgbXlfdGhlbWUoMTIpICsNCiAgICAgdGhlbWUoYXhpcy50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSAxMiwgY29sb3IgPSAiYmxhY2siKSwNCiAgICAgYXhpcy50aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gMTUsIGNvbG9yID0gImJsYWNrIikpIA0KZy50cnQgDQpgYGANCg0KQ3Xhu5FpIGPDuW5nLCB0YSBjw7MgdGjhu4Mga2jhuqNvIHPDoXQga2h1eW5oIGjGsOG7m25nIMSRw6FwIOG7qW5nIMSRaeG7gXUgdHLhu4sg4bufIG3hu5dpIGPDoSB0aOG7gyBi4buHbmggbmjDom4sIGLDonkgZ2nhu50gdGjDrCBi4bqhbiBjw7MgdGjhu4MgaGnhu4N1IHThuqFpIHNhb3RhIGTDuW5nIG3DtCBow6xuaCB24bubaSByYW5kb20gZWZmZWN0DQoNCmBgYHtyLG1lc3NhZ2UgPSBGQUxTRSx3YXJuaW5nPUZBTFNFfQ0KZGF0JT4lZ2dwbG90KGFlcyh4PXBlcmlvZCx5PXNlaSxncm91cD0xKSkrDQogIGdlb21fcGF0aChhZXMoY29sb3I9dHJ0KSxhbHBoYT0wLjksc2l6ZT0xKSsNCiAgZ2VvbV9wb2ludChhZXMoY29sb3I9dHJ0LHNoYXBlPXRydCksYWxwaGE9MC45LHNpemU9MikrDQogIG15X3RoZW1lKDUpKw0KICBmYWNldF93cmFwKH5zdWJqZWN0LHNjYWxlcz0iZnJlZSIsbmNvbD0xMCkNCg0KYGBgDQoNCiMgTcO0IGjDrG5oIE5CSUkgQmF5ZXMgY2jDrW5oIHRo4bupYw0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NClhtYXQgPC0gIG1vZGVsLm1hdHJpeCh+ZGF0JHRydCtkYXQkbGJhc2UrZGF0JHZpc2l0KQ0KTnN1YiA8LSAgbmxldmVscyhkYXQkc3ViamVjdCkgICMgbm8gb2Ygc3ViamVjdHMNCksgPSBuY29sKFhtYXQpDQpkYXRhLk5CSSA9IGxpc3QoeSA9IGRhdCRzZWksICAgICAgICAgICAjICByZXNwb25zZSB2YXJpYWJsZQ0KICAgICAgICAgICAgICAgIE4gPSBucm93KGRhdCksICAgICAgICAgICMgbm8gb2Ygb2JzDQogICAgICAgICAgICAgICAgWCA9IFhtYXQsICAgICAgICAgICAgICAgIyBtb2RlbCBtYXRyaXggb2YgcHJlZGljdG9ycyBhbmQgaW50ZXJjZXB0DQogICAgICAgICAgICAgICAgSyA9IEssICAgICAgICAgICMgbm8gb2YgcHJlZGljdG9ycyArIGludGVyY2VwdA0KICAgICAgICAgICAgICAgIE5zdWIgPSBOc3ViLCAgICAgICAgICAgICMgbm8gb2Ygc3VicyAoIDU5ID0gbm8gb2YgcmFuZG9tIGVmZmVjdCBWQUxVRXMpDQogICAgICAgICAgICAgICAgc3ViID0gZGF0JHN1YmplY3QgJT4lIGFzLm51bWVyaWMoLiksICMgU3ViamVjdCBJRA0KICAgICAgICAgICAgICAgIGJwcmkgPSByZXAoMCwgSyksICAgICAgICMgVG8gZGVmaW5lIHByaW9yIGZvciBiZXRhDQogICAgICAgICAgICAgICAgQnByaSA9IGRpYWcoMSwgSyksICAgICAgIyBQcmlvciBtYXRyaXggb2YgYmV0YSANCiAgICAgICAgICAgICAgICBhcHJpID0gcmVwKDAsIE5zdWIpLCAgICAjIFRvIGRlZmluZSBwcmlvciBmb3IgaW50ZXJjZXB0DQogICAgICAgICAgICAgICAgQXByaSA9IGRpYWcoMSwgTnN1YiApKSAgIyBQcmlvciBtYXRyaXggb2YgaW50ZXJjZXB0IA0KDQpgYGANCg0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmZpdC5zdGFuIDwtICBzdGFuKG1vZGVsX2NvZGUgPSBzdGFuX05CSS5zdHJpLCANCiAgICAgICAgICAgICAgICAgIGRhdGEgPSBkYXRhLk5CSSwgDQogICAgICAgICAgICAgICAgICBjaGFpbnM9NCwgDQogICAgICAgICAgICAgICAgICBpdGVyPTMwMDAsIA0KICAgICAgICAgICAgICAgICAgd2FybXVwPTEwMDAsIA0KICAgICAgICAgICAgICAgICAgdGhpbj0yKQ0KYGBgDQoNCkvhur90IHF14bqjIGPhu6dhIG3DtCBow6xuaCBCYXllcyBOQkkyIMSRxrDhu6NjIHTDs20gdOG6r3QgbmjGsCBzYXU6DQoNCmBgYHtyfQ0KDQpwcmludChmaXQuc3RhbiwgZGlnaXRzPTMsIHBhcnM9YygiYmV0YSIsInNpZ21hX3JlIiwiYWxwaGEiKSkNCg0KYGBgDQoNClRhIGNo4buJIHF1YW4gdMOibSDEkeG6v24gcGjDom4gcGjhu5FpIGjhuq11IMSR4buLbmggY+G7p2EgMyBoaeG7h3Ug4bupbmcgY2jDrW5oLGTGsOG7m2kgZOG6oW5nIEluY2lkZW5jZSByYXRlIHJhdGlvIChJUlIgPSBleHAoYmV0YSkpDQoNCmBgYHtyLG1lc3NhZ2UgPSBGQUxTRSx3YXJuaW5nPUZBTFNFfQ0KcG9zdGRmPWFzLmRhdGEuZnJhbWUoZml0LnN0YW4pJT4lYXNfdGliYmxlKCkNCg0KYmV0YWRmPXBvc3RkZlssYygxOjQpXQ0KDQpuYW1lcyhiZXRhZGYpPWMoIkludGVyY2VwdCIsIlByb2dhYmlkZSIsIkxvZ0Jhc2UiLCJWaXNpdCIpDQoNCmJldGFkZj1iZXRhZGYlPiVtdXRhdGUoLixJdGVyYXRpb249YXMubnVtZXJpYyhyZXAoYygxOjEwMDApLDQpKSwNCiAgICAgICAgICAgICAgICAgICAgICAgIENoYWluPWFzLmZhY3RvcihyZXAoYygxOjQpLGVhY2g9MTAwMCkpKQ0KDQpoZWFkKGJldGFkZiklPiUuWyxjKDE6NCldJT4la25pdHI6OmthYmxlKCkNCg0KcDE9YmV0YWRmJT4lZ2F0aGVyKEludGVyY2VwdDpWaXNpdCxrZXk9IlN0ZXAiLHZhbHVlPSJmaXgiKSU+JQ0KICBnZ3Bsb3QoYWVzKHg9ZXhwKGZpeCksZmlsbD1DaGFpbixjb2w9Q2hhaW4pKSsNCiAgZ2VvbV9kZW5zaXR5KGFscGhhPTAuMyxzaG93LmxlZ2VuZCA9IEYpKw0KICBmYWNldF93cmFwKH5TdGVwLG5jb2w9MSxzY2FsZXMgPSAiZnJlZV95IikrDQogIHNjYWxlX3hfY29udGludW91cygiSVJSIikrDQogIHRoZW1lX2J3KDYpDQoNCnAyPWJldGFkZiU+JWdhdGhlcihJbnRlcmNlcHQ6VmlzaXQsa2V5PSJTdGVwIix2YWx1ZT0iZml4IiklPiUNCiAgZ2dwbG90KGFlcyh5PWV4cChmaXgpLHg9SXRlcmF0aW9uLGNvbD1DaGFpbikpKw0KICBnZW9tX3BhdGgoYWxwaGE9MC41LHNob3cubGVnZW5kID0gRikrDQogIGZhY2V0X3dyYXAoflN0ZXAsbmNvbD0xLHNjYWxlcyA9ICJmcmVlIikrDQogIHNjYWxlX3lfY29udGludW91cygiSVJSIikrDQogIHRoZW1lX2J3KDYpDQoNCmdyaWRFeHRyYTo6Z3JpZC5hcnJhbmdlKHAyLHAxLG5jb2w9MikNCmBgYA0KDQpO4bq/dSBjaOG7iSB4w6l0IEluY2lkZW5jZSByYXRlIGNoYW5nZSBjaG8geeG6v3UgdOG7kSDEkWnhu4F1IHRy4buLLCB0YSBjw7MgdGjhu4MgZGnhu4VuIGdp4bqjaTpUaHXhu5FjIHByb2dhYmlkIGzDoG0gZ2nhuqNtIHThuqduIHN14bqldCBjxqFuIMSR4buZbmcga2luaCBraG/huqNuZyAyNi43MiAlICg0LjM2ICUgxJHhur9uIDQzLjclICkgc28gduG7m2kgbmjDs20gUGxhY2Viby4NCg0KYGBge3IsbWVzc2FnZSA9IEZBTFNFLHdhcm5pbmc9RkFMU0V9DQpIbWlzYzo6ZGVzY3JpYmUoZXhwKGJldGFkZiRQcm9nYWJpZGUpLTEpDQpgYGANCg0KVMawxqFuZyBxdWFuIGdp4buvYSBJbmNpZGVuY2UgcmF0ZSBjaGFuZ2UgY+G7p2EgVmlzaXQgdsOgIFByb2dhYmlkZSDEkcaw4bujYyBiaeG7g3UgZGnhu4VuIG5oxrAgc2F1Og0KDQpgYGB7cixtZXNzYWdlID0gRkFMU0Usd2FybmluZz1GQUxTRX0NCmJldGFkZiRwc2V1ZG9Hcm91cD1mYWN0b3IocmVwKGMoMToyMCksZT1ucm93KGJldGFkZikvMjApKQ0KDQpiZXRhZGYlPiVnYXRoZXIoVmlzaXQsUHJvZ2FiaWRlLGtleT0iRmFjdG9ycyIsdmFsdWU9ImZpeCIpJT4lDQogIGdncGxvdChhZXMoeT1GYWN0b3JzLA0KICAgICAgICAgICAgIHg9MTAwKihleHAoZml4KS0xKQ0KICAgICAgICAgICAgICkpKw0KICBnZW9tX2RlbnNpdHlfcmlkZ2VzX2dyYWRpZW50KGFlcyhmaWxsID0gLi54Li4pLA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHNob3cubGVnZW5kID0gRiwNCiAgICAgICAgICAgc2NhbGU9MC44LA0KICAgICAgICAgICBncmFkaWVudF9sd2QgPSAwLjUpKw0KICBnZW9tX2RlbnNpdHlfcmlkZ2VzKGFlcyhjb2w9cHNldWRvR3JvdXApLA0KICAgICAgICAgICAgICAgICAgICAgIHNjYWxlPTAuOSxhbHBoYT0wLjAxLA0KICAgICAgICAgICAgICAgICAgICAgIHNob3cubGVnZW5kID0gRikrDQogIGdlb21fdmxpbmUoeGludGVyY2VwdD1jKDAsLTEwLDEwKSxjb2w9YygicmVkMyIsImJsdWUzIiwiYmxhY2siKSxsaW5ldHlwZT0yKSsNCiAgc2NhbGVfeF9jb250aW51b3VzKCJJbmNpZGVuY2UgcmF0ZSBjaGFuZ2UgKCUpIiwNCiAgICAgICAgICAgICAgICAgICAgIGV4cGFuZCA9IGMoMC4wMSwgMCksDQogICAgICAgICAgICAgICAgICAgICBicmVha3M9YygtNzAsLTYwLC01MCwtNDAsLTMwLC0yMCwtMTAsMCwxMCwyMCwzMCw0MCw1MCkpKw0KICBzY2FsZV95X2Rpc2NyZXRlKGV4cGFuZCA9IGMoMC4wMSwgMCkpKw0KICBzY2FsZV9maWxsX3ZpcmlkaXMob3B0aW9uID0gIkEiLGJlZ2luPTAuMixlbmQ9MSkrDQogIHNjYWxlX2NvbG91cl92aXJpZGlzKG9wdGlvbiA9ICJBIixkaXNjcmV0ZSA9IFQpKw0KICB0aGVtZV9idygpDQogIA0KYGBgDQoNClBow6JuIHTDrWNoIFJPUEUgdGhlbyBKLkpydXNjaGtlIHbhu5tpIG5nxrDhu6FuZyB0csOqbiB2w6AgZMaw4bubaSBs4bqnbiBsxrDhu6N0IGzDoCAtNSUgdsOgICs1JSAsIGPDuW5nIHBow6JuIHTDrWNoIENvbXBWYWwgduG7m2kgbmfGsOG7oW5nIGzDoCAxMCUNCg0KYGBge3IsbWVzc2FnZSA9IEZBTFNFLHdhcm5pbmc9RkFMU0V9DQpIRElGPSBmdW5jdGlvbiggc2FtcGxlVmVjLGNyZWRNYXNzPTAuOTc1ICkgew0KICBzb3J0ZWRQdHMgPSBzb3J0KCBzYW1wbGVWZWMgKQ0KICBjaUlkeEluYyA9IGNlaWxpbmcoIGNyZWRNYXNzICogbGVuZ3RoKCBzb3J0ZWRQdHMgKSApDQogIG5DSXMgPSBsZW5ndGgoIHNvcnRlZFB0cyApIC0gY2lJZHhJbmMNCiAgY2lXaWR0aCA9IHJlcCggMCAsIG5DSXMgKQ0KICBmb3IgKCBpIGluIDE6bkNJcyApIHsNCiAgICBjaVdpZHRoWyBpIF0gPSBzb3J0ZWRQdHNbIGkgKyBjaUlkeEluYyBdIC0gc29ydGVkUHRzWyBpIF0NCiAgfQ0KICBIREltaW4gPSBzb3J0ZWRQdHNbIHdoaWNoLm1pbiggY2lXaWR0aCApIF0NCiAgSERJbWF4ID0gc29ydGVkUHRzWyB3aGljaC5taW4oIGNpV2lkdGggKSArIGNpSWR4SW5jIF0NCiAgSERJbGltID0gYyggSERJbWluICwgSERJbWF4ICkNCiAgcmV0dXJuKCBIRElsaW0gKQ0KfQ0KDQpTVU1LPWZ1bmN0aW9uKHBhcmFtU2FtcGxlVmVjLGNvbXBWYWw9TlVMTCAsIFJPUEU9TlVMTCAsIGNyZWRNYXNzPTAuOTc1KSB7DQogIG1lYW5QYXJhbSA9IG1lYW4oIHBhcmFtU2FtcGxlVmVjICkNCiAgbWVkaWFuUGFyYW0gPSBtZWRpYW4oIHBhcmFtU2FtcGxlVmVjICkNCiAgZHJlcyA9IGRlbnNpdHkoIHBhcmFtU2FtcGxlVmVjICkNCiAgbW9kZVBhcmFtID0gZHJlcyR4W3doaWNoLm1heChkcmVzJHkpXQ0KICBoZGlMaW0gPSBIRElGKCBwYXJhbVNhbXBsZVZlYyAsIGNyZWRNYXNzPWNyZWRNYXNzICkNCiAgaWYgKCAhaXMubnVsbChjb21wVmFsKSApIHsNCiAgICBwY2d0Q29tcFZhbCA9ICggMTAwICogc3VtKCBwYXJhbVNhbXBsZVZlYyA+IGNvbXBWYWwgKSANCiAgICAgICAgICAgICAgICAgICAgLyBsZW5ndGgoIHBhcmFtU2FtcGxlVmVjICkgKQ0KICB9IGVsc2Ugew0KICAgIGNvbXBWYWw9TkENCiAgICBwY2d0Q29tcFZhbD1OQQ0KICB9DQogIGlmICggIWlzLm51bGwoUk9QRSkgKSB7DQogICAgcGNsdFJvcGUgPSAoIDEwMCAqIHN1bSggcGFyYW1TYW1wbGVWZWMgPCBST1BFWzFdICkgDQogICAgICAgICAgICAgICAgIC8gbGVuZ3RoKCBwYXJhbVNhbXBsZVZlYyApICkNCiAgICBwY2d0Um9wZSA9ICggMTAwICogc3VtKCBwYXJhbVNhbXBsZVZlYyA+IFJPUEVbMl0gKSANCiAgICAgICAgICAgICAgICAgLyBsZW5ndGgoIHBhcmFtU2FtcGxlVmVjICkgKQ0KICAgIHBjaW5Sb3BlID0gMTAwLShwY2x0Um9wZStwY2d0Um9wZSkNCiAgfSBlbHNlIHsgDQogICAgUk9QRSA9IGMoTkEsTkEpDQogICAgcGNsdFJvcGU9TkEgDQogICAgcGNndFJvcGU9TkEgDQogICAgcGNpblJvcGU9TkEgDQogIH0gIA0KICByZXR1cm4oIGMoIE1lYW49bWVhblBhcmFtICwgTWVkaWFuPW1lZGlhblBhcmFtICwgTW9kZT1tb2RlUGFyYW0gLCANCiAgICAgICAgICAgICBIRElsZXZlbD1jcmVkTWFzcyAsIExMPWhkaUxpbVsxXSAsIFVMPWhkaUxpbVsyXSAsIA0KICAgICAgICAgICAgIENvbXBWYWw9Y29tcFZhbCAsIFBjbnRHdENvbXBWYWw9cGNndENvbXBWYWwgLCANCiAgICAgICAgICAgICBST1BFbG93PVJPUEVbMV0gLCBST1BFaGlnaD1ST1BFWzJdICwNCiAgICAgICAgICAgICBQY250THRST1BFPXBjbHRSb3BlICwgUGNudEluUk9QRT1wY2luUm9wZSAsIFBjbnRHdFJPUEU9cGNndFJvcGUgKSApDQp9DQoNCnN1bW1hcnlLcnVzY2hrZT1mdW5jdGlvbihNQ01DLGNvbXBWYWw9TlVMTCwgcm9wZT1OVUxMLGNyZWRNYXNzPU5VTEwpew0KICBzdW1tYXJ5SW5mbyA9IE5VTEwNCiAgc3VtbWFyeUluZm8gPSBjYmluZChzdW1tYXJ5SW5mbywgIkVzdGltYXRlZCI9IFNVTUsoTUNNQywNCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY29tcFZhbD1jb21wVmFsLA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBST1BFPXJvcGUsY3JlZE1hc3M9Y3JlZE1hc3MpKQ0KICByZXR1cm4oc3VtbWFyeUluZm8pDQp9DQoNCklSQz0xMDAqKDEtZXhwKGJldGFkZiRQcm9nYWJpZGUpKSU+JWFzX2RhdGFfZnJhbWUoKQ0KDQpzdW1tYXJ5S3J1c2Noa2UoTUNNQyA9IElSQyR2YWx1ZSwNCiAgICAgICAgICAgICAgICAgICAgICAgY29tcFZhbD0xMCwNCiAgICAgICAgICAgICAgICAgICAgICAgcm9wZT1jKC01LDUpLA0KICAgICAgICAgICAgICAgICAgICAgICBjcmVkTWFzcz0wLjk3NSklPiVhcy5kYXRhLmZyYW1lKCkNCg0KYmV0YWRmJElSQz0xMDAqKDEtZXhwKGJldGFkZiRQcm9nYWJpZGUpKQ0KYGBgDQogDQogS+G6v3QgcXXhuqMgY2hvIHRo4bqleTogOTQuNiUgbeG6rXQgxJHhu5kgcGjDom4gcGjhu5FpIGjhuq11IMSR4buLbmggbuG6sW0gbmdvw6BpIG5nxrDhu6FuZyB0aGF5IMSR4buVaSA1JSwgY2jhu4kgY8OzIDQuMDc1JSBt4bqtdCDEkeG7mSByxqFpIHbDoG8gdHJvbmcga2hv4bqjbmcgKy8tNSU7IDg5Ljk1JSBt4bqtdCDEkeG7mSBu4bqxbSBuZ2/DoGkgbmfGsOG7oW5nIDEwJQ0KIA0KxJDhu4MgY2jhuq9jIGNo4bqvbiwgdGEgbMOgbSB0aMOqbSAxIHBow6JuIHTDrWNoIEJheWVzIEZhY3RvciBjaG8gMjEgbmfGsOG7oW5nICB0aGF5IMSR4buVaSBj4bunYSBJUiB04burIDAgxJHDqW4tMjAlDQoNCmBgYHtyfQ0KDQpsaWJyYXJ5KGJybXMpDQoNCmJldGFkZiRJUkM9MTAwKigxLWV4cChiZXRhZGYkUHJvZ2FiaWRlKSkNCg0KdGhyZXNob2xkPXJlcChOQSwyMSkNCkJheWVzRmFjdG9yPXJlcChOQSwyMSkNCg0KdGhyZXM9YygwOjIwKQ0KDQpmb3IoaSBpbiAoMToyMSkpew0KICB0aHI9dGhyZXNbaV0NCiAgdGhyZXNob2xkW2ldPXRocg0KICBoeXA9cGFzdGUoIklSQz4iLHRocixzZXA9IiIpDQogIGJmPWJybXM6Omh5cG90aGVzaXMoYmV0YWRmLGh5cCxhbHBoYT0wLjA1KQ0KICBCYXllc0ZhY3RvcltpXT1iZiRoeXBvdGhlc2lzJEV2aWQuUmF0aW8NCn0NCg0KYmZkZj1jYmluZCh0aHJlc2hvbGQsQmF5ZXNGYWN0b3IpJT4lYXNfdGliYmxlKCkNCg0KYmZkZiU+JWdncGxvdChhZXMoeD10aHJlc2hvbGQseT1CYXllc0ZhY3RvcixmaWxsPUJheWVzRmFjdG9yKSkrDQogIGdlb21fcGF0aCgpKw0KICBnZW9tX3BvaW50KHNob3cubGVnZW5kID0gRixzaXplPTUsc2hhcGU9MjEsY29sPSJibGFjayIpKw0KICBnZW9tX3RleHQoYWVzKGxhYmVsPXJvdW5kKEJheWVzRmFjdG9yLDIpKSxjb2w9ImJsYWNrIixzaG93LmxlZ2VuZCA9IEYsYW5nbGUgPSA2MCxudWRnZV95PXNxcnQoYmZkZiRCYXllc0ZhY3RvcikrMSxudWRnZV94PTAsc2l6ZT00KSsNCiAgdGhlbWVfYncoKStzY2FsZV94X2NvbnRpbnVvdXMoYnJlYWtzPWMoMDoyMCkpKw0KICBnZW9tX2hsaW5lKHlpbnRlcmNlcHQgPSBjKDEwLDMwKSxsaW5ldHlwZT0yLGNvbD0iYmx1ZSIpKw0KICBzY2FsZV9maWxsX3ZpcmlkaXMob3B0aW9uPSJEIixkaXJlY3Rpb249LTEpDQpgYGANCg0KTuG6v3UgbOG6pXkgbmfGsOG7oW5nIGto4bqjIHTDrW4gY+G7p2EgQkY9MzAsIHRhIGPDsyB0aOG7gyBr4bq/dCBsdeG6rW4gbMOgIHR1eSBQcm9nYWJpZGUgY8OzIGPhuqNpIHRoaeG7h24gdHJp4buHdSBjaOG7qW5nIMSR4buZbmcga2luaCwgbmjGsG5nIGhp4buHdSDhu6luZyBuw6B5IGtow6Egbmjhu48sIFbhu5tpIEgxIGzDoCBJUiB0aGF5IMSR4buVaSA+IDEwJSB0aMOsIEJGIGNo4buJIGPDsyA4Ljk1OyBUYSBjaOG7iSBjw7MgdGjhu4MgdGluIGPhuq15IGhp4buHdSDhu6luZyDEkWnhu4F1IHRy4buLIHThu6sgMS01JSBtw6AgdGjDtGkuDQogDQojIERp4buFbiDEkeG6oXQgdsSDbiBi4bqjbiBraG9hIGjhu41jDQoNClBoxrDGoW5nIHBow6FwIHRo4buRbmcga8OqIDoNCg0KRG8ga+G6v3QgcXXhuqMgKHPhu5EgY8ahbiDEkeG7mW5nIGtpbmgpIGzDoCBt4buZdCBiaeG6v24gc+G7kSDEkeG6v20gdsOgIGLhu4sgcGjDom4gdMOhbiwgbsOzIMSRxrDhu6NjIGto4bqjbyBzw6F0IHRoZW8gcXV5IGx14bqtdCBwaMOibiBwaOG7kWkgbmjhu4sgdGjhu6ljIMOibSB0eXBlIElJIChOQklJKS4gSGnhu4d1IOG7qW5nIGPhu6dhIHRodeG7kWMgUHJvZ2FiaWRlIHNvIHbhu5tpIFBsYWNlYm8gbMOqbiBr4bq/dCBj4bulYyBsw6JtIHPDoG5nIMSRxrDhu6NjIGto4bqjbyBzw6F0IGLhurFuZyBt4buZdCBtw7QgaMOsbmggaOG7k2kgcXV5IGjhu5duIGjhu6NwIChNaXhlZCBtb2RlbCkgdGhlbyBwaMawxqFuZyBwaMOhcCBCYXllcy4gTcO0IGjDrG5oIG7DoHkgY8OzIGhp4buHdSBjaOG7iW5oIGNobyB0w6xuaCB0cuG6oW5nIGLhu4duaCB04bqhaSBraOG7n2kgxJFp4buDbSBuZ2hpw6puIGPhu6l1IChsYmFzZSksIHnhur91IHThu5EgdGjhu51pIGdpYW4gKDQgdHLhu41uZyBz4buRIHTGsMahbmcgcGjhuqNuKSB2w6AgaGnhu4d1IOG7qW5nIG5n4bqrdSBuaGnDqm4gY8OhIHRo4buDLiBUw61uaCBo4bujcCBsw70gY+G7p2EgbcO0IGjDrG5oIMSRxrDhu6NjIMSRw6FuaCBnacOhIGLhurFuZyB0acOqdSBjaMOtIFdBSUMuIFN1eSBkaeG7hW4gdGjhu5FuZyBrw6ogZOG7sWEgdsOgbyBt4bqtdCDEkeG7mSBwaMOibiBwaOG7kWkgaOG6rXUgxJHhu4tuaCB2w6AgdOG7iSB0cuG7jW5nIGNo4bupbmcgY+G7qSAoQmF5ZXMgZmFjdG9yKSB24bubaSBnaeG6oyB0aHV54bq/dCBIMSBsw6AgdGh14buRYyBQcm9nYWJpZGUgZ2nhuqNtIHThuqduIHN14bqldCDEkeG7mW5nIGtpbmggw610IG5o4bqldCBsw6AgNSUgc28gduG7m2kgUGxhY2Viby4NCg0KS+G6v3QgcXXhuqM6DQoNClBow6JuIHBo4buRaSBo4bqtdSDEkeG7i25oIGPhu6dhIHPhu7EgdGhheSDEkeG7lWkgdOG6p24gc3XhuqV0IHBow6F0IHNpbmggY8ahbiDEkeG7mW5nIGtpbmggKElSIGNoYW5nZSwlKSBjaG8gdGjhuqV5IHRodeG7kWMgUHJvZ2FiaWQgY+G6o2kgdGhp4buHbiB0cnVuZyBiw6xuaCAyNi43MiUgc+G7kSBjxqFuIMSR4buZbmcga2luaCBzbyB24bubaSBuaMOzbSBQbGFjZWJvLiDEkOG6v24gOTQuNiUgbeG6rXQgxJHhu5kgcGjDom4gcGjhu5FpIGjhuq11IMSR4buLbmggbuG6sW0gbmdvw6BpIG5nxrDhu6FuZyB0aGF5IMSR4buVaSA1JSwgY2jhu4kgY8OzIDQuMDc1JSBt4bqtdCDEkeG7mSByxqFpIHbDoG8gdHJvbmcga2hv4bqjbmcgKy8tNSU7IDg5Ljk1JSBt4bqtdCDEkeG7mSBu4bqxbSBuZ2/DoGkgbmfGsOG7oW5nIDEwLiBUdXkgbmhpw6puIEJheWVzIEZhY3RvciBjaG8gdGjhuqV5IGto4bqjIG7Eg25nIHjDoWMgdMOtbiB24buBIGhp4buHdSBxdeG6oyBuw6B5IGzDoCBraMOhIHRo4bqlcCAoQkY9MTcuNTIgY2hvIG5nxrDhu6FuZyA1JSB2w6AgOC45NSBjaG8gbmfGsOG7oW5nIDEwJSkuDQoNCiMgVOG7lW5nIGvhur90DQoNCkLDoGkgdGjhu7FjIGjDoG5oIHPhu5EgNyDEkeG6v24gxJHDonkgbMOgIGNo4bqlbSBk4bupdC4gVHJvbmcgYsOgaSBuw6B5IGNow7puZyB0w7RpIMSRw6MgdHJ1eeG7gW4gdOG6o2kgbeG7mXQgc+G7kSB0aMO0bmcgxJFp4buHcCBjaMOtbmggbmjGsCBzYXU6DQoNCjEpIEtow6FjIHbhu5tpIHRow60gbmdoaeG7h20geSBo4buNYyBsw6JtIHPDoG5nLCDEkeG7kWkgdMaw4bujbmcgY+G7p2EgcGjDom4gdMOtY2ggdGjhu5FuZyBrw6oga2jDtG5nIHBo4bqjaSBsw6Agbmjhu69uZyDEkeG6oWkgbMaw4bujbmcgc2luaCBsw70gYuG7h25oLCBtw6AgbMOgIG5o4buvbmcgYmnhur9uIHPhu5Egbmfhuqt1IG5oacOqbiDEkeG6oWkgZGnhu4duIGNobyBjw6FjIMSR4bqhaSBsxrDhu6NuZyBuw6B5LiBWaeG7h2MgbOG7sWEgY2jhu41uIHF1eSBsdeG6rXQgcGjDom4gcGjhu5FpIHBow7kgaOG7o3AgxJHhu4MgY2h1eeG7g24gdOG7qyDEkeG6oWkgbMaw4bujbmcgc2FuZyBiaeG6v24gc+G7kSB5w6p1IGPhuqd1IG5nxrDhu51pIGLDoWMgc8SpIHBo4bqjaSB24bqtbiBk4bulbmcgbMO9IHRodXnhur90IHjDoWMgc3XhuqV0Lg0KDQoyKSBOaOG7r25nIGJp4bq/biBz4buRIMSR4bq/bSBraMO0bmcgdGjhu4MgxJHGsOG7o2Mga2jhuqNvIHPDoXQgYuG6sW5nIHBoxrDGoW5nIHBow6FwIHRydXnhu4FuIHRo4buRbmcgZOG7sWEgdHLDqm4gZ2nhuqMgxJHhu4tuaCBwaMOibiBwaOG7kWkgR2F1c3NpYW4uIE3hu5l0IHPhu5EgZ2nhuqNpIHBow6FwIMSRw6MgxJHGsOG7o2MgZ2nhu5tpIHRoaeG7h3UsIGJhbyBn4buTbSBtw7QgaMOsbmggUG9pc3NvbiwgUXVhc2ktUG9pc3NvbiB2w6AgTmjhu4sgdGjhu6ljIMOibSB0eXBlIDIuDQoNCjMpIFBow6JuIHTDrWNoIEJheWVzIHRyb25nIFNUQU4gY2hvIHBow6lwIGThu7FuZyBtw7QgaMOsbmggR0xNIHbhu5tpIHBow6JuIHBo4buRaSB1eeG7g24gY2h1eeG7g24sIG5oxrAgdHJvbmcgdHLGsOG7nW5nIGjhu6NwIG7DoHkgbMOgIHBow6JuIHBo4buRaSBOQklJLiBTdXkgZGnhu4VuIEJheWVzIGPFqW5nIG1hbmcgbOG6oWkgbmhp4buBdSB0aMO0bmcgdGluIGjGoW4gY8OhYyBraeG7g20gxJHhu4tuaCBj4buVIMSRaeG7g24uDQoNCk5ow7NtIEJBViB4aW4gY2jDom4gdGjDoG5oIGPhuqNtIMahbiBz4buxIOG7p25nIGjhu5kgY+G7p2EgY8OhYyBi4bqhbi4gQ2jDum5nIHTDtGkgc+G6vSBn4bq3cCBs4bqhaSBjw6FjIGLhuqFuIHRyb25nIG3hu5l0IGLDoGkga2jDoWMuIA==