From c3b41446898fa328c71b6d7f9a1ca49c1d4b2552 Mon Sep 17 00:00:00 2001 From: Louis Lacoste Date: Mon, 22 Apr 2024 10:26:30 +0200 Subject: [PATCH] Changing analysis --- code/analysis/analyze_NA_robustness.R | 96 --------- code/analysis/analyze_NA_robustness.Rmd | 183 ++++++++++++++++++ ...-04-2024_15-18-55_uniform_nested_1-200.Rds | Bin 12553 -> 0 bytes 3 files changed, 183 insertions(+), 96 deletions(-) delete mode 100644 code/analysis/analyze_NA_robustness.R create mode 100644 code/analysis/analyze_NA_robustness.Rmd delete mode 100644 code/results/simulations/NA_robustness/NA_robustness_19-04-2024_15-18-55_uniform_nested_1-200.Rds diff --git a/code/analysis/analyze_NA_robustness.R b/code/analysis/analyze_NA_robustness.R deleted file mode 100644 index f01e7fb..0000000 --- a/code/analysis/analyze_NA_robustness.R +++ /dev/null @@ -1,96 +0,0 @@ -library(ggplot2) -library(tidyverse) -library(latex2exp) -library(patchwork) - -#  Loading data -data_folder <- file.path("code", "results", "simulations", "NA_robustness") - -data <- readRDS(file.path( - data_folder, - "NA_robustness_19-04-2024_15-18-55_uniform_nested_1-200.Rds" -)) - -data[["model"]] <- factor(data[["model"]], levels = c("iid", "pi", "rho", "pirho")) - -max_repetition <- max(data$repetition) - -#  Averaging over repetitions -averaged_data <- data %>% - group_by(prop_NAs, model) %>% - select(-repetition) %>% - summarise_all(list(mean = mean, sd = sd)) - -#  Preparing auc_data -auc_data <- averaged_data %>% - select(c(prop_NAs, model) | contains("auc_")) %>% - rename_with(~ gsub("auc_", "", .x, fixed = TRUE)) %>% - filter(prop_NAs != 0) - -auc_data_long <- - bind_cols( - auc_data %>% select(c("prop_NAs", "model") | contains("_mean")) %>% - pivot_longer(!c(prop_NAs, model), - names_to = c("method"), - values_to = "auc_mean" - ), - auc_data %>% select(c("prop_NAs", "model") | contains("_sd")) %>% - pivot_longer(!c(prop_NAs, model), - names_to = NULL, - values_to = "auc_sd" - ) %>% ungroup() %>% select(!c("prop_NAs", "model")) - ) %>% mutate(method = method %>% gsub( - pattern = "_mean", - replacement = "", fixed = TRUE - )) -# Preparing ARI data -ari_data <- averaged_data %>% select(c(prop_NAs, model) | contains("ari")) - -ari_data_long <- bind_cols( - ari_data %>% select(c("prop_NAs", "model") | contains("_mean")) %>% - pivot_longer(!c(prop_NAs, model), - names_to = c("method"), - values_to = "ari_mean" - ), - ari_data %>% select(c("prop_NAs", "model") | contains("_sd")) %>% - pivot_longer(!c(prop_NAs, model), - names_to = NULL, - values_to = "ari_sd" - ) %>% ungroup() %>% select(!c("prop_NAs", "model")) -) %>% - separate(method, into = c("dim", "method", "useless"), sep = "_") %>% - select(-useless) - -auc_plot <- ggplot(auc_data_long) + - aes(x = prop_NAs, y = auc_mean) + - geom_line(aes(color = method)) + - geom_point(aes(color = method)) + - geom_ribbon(aes(ymin = auc_mean - auc_sd, ymax = auc_mean + auc_sd, fill = method), alpha = 0.2) + - ylim(c(min(auc_data_long[["auc_mean"]]), max(auc_data_long[["auc_mean"]]))) + - scale_x_continuous(breaks = scales::pretty_breaks(n = 10L)) + - ylab(TeX("\\bar{AUC}")) + - xlab("NA proportion") + - labs(fill = "Method", color = "Method") + - ggtitle(TeX(paste( - "$\\bar{AUC}\\pm s_{AUC}$", ", function of NA proportion. N=", max_repetition - ))) + - facet_grid(rows = vars(model), as.table = TRUE) + - theme_bw() - - -ari_plot <- ggplot(ari_data_long) + - aes(x = prop_NAs, y = ari_mean) + - geom_line(aes(color = method)) + - geom_point(aes(color = method)) + - geom_ribbon(aes(ymin = ari_mean - ari_sd, ymax = ari_mean + ari_sd, fill = method), alpha = 0.2) + - scale_x_continuous(breaks = scales::pretty_breaks(n = 10L)) + - ylab(TeX("$\\bar{ARI^d}$")) + - xlab("NA proportion") + - labs(fill = "Method", color = "Method") + - ggtitle(TeX(paste( - "$\\bar{ARI^d}\\pm s_{ARI^d}$", ", function of NA proportion. N=", max_repetition - ))) + - facet_grid(rows = vars(model), cols = vars(dim)) + - theme_bw() - -(auc_plot | ari_plot) + patchwork::plot_layout(guides = "collect") diff --git a/code/analysis/analyze_NA_robustness.Rmd b/code/analysis/analyze_NA_robustness.Rmd new file mode 100644 index 0000000..cc40a60 --- /dev/null +++ b/code/analysis/analyze_NA_robustness.Rmd @@ -0,0 +1,183 @@ +--- +title: "Analysis NA robustness" +output: + bookdown::html_document2: + toc: true + fig_caption: true +--- + +```{r knitropts, echo = FALSE} +knitr::opts_chunk$set(fig.width=12L, fig.height=10L) +``` + +```{r imports, echo = FALSE, include=FALSE} +suppressPackageStartupMessages(library(ggplot2)) +suppressPackageStartupMessages(library(tidyverse)) +suppressPackageStartupMessages(library(latex2exp)) +suppressPackageStartupMessages(library(here)) +loadNamespace(package = "patchwork") +``` + +```{r data, echo = FALSE} +#  Loading data +data_folder <- file.path(here(), "code", "results", "simulations", "NA_robustness") +files_list <- list.files(data_folder, pattern = "^NA_robustness_19-04-2024") +``` +```{r df_func, echo = FALSE} +prepare_dataframes_auc_ari <- function(data) { + data[["model"]] <- factor(data[["model"]], levels = c("iid", "pi", "rho", "pirho")) + + max_repetition <- max(data[["repetition"]]) + sampling <- unique(data[["sampling"]]) + struct <- unique(data[["struct"]]) + + #  Averaging over repetitions + averaged_data <- data %>% + group_by(prop_NAs, model) %>% + select(-c(repetition, sampling, struct)) %>% + summarise_all(list(mean = mean, sd = sd)) + + #  Preparing auc_data + auc_data <- averaged_data %>% + select(c(prop_NAs, model) | contains("auc_")) %>% + rename_with(~ gsub("auc_", "", .x, fixed = TRUE)) %>% + filter(prop_NAs != 0) + + auc_data_long <- + bind_cols( + auc_data %>% + select(c("prop_NAs", "model") | contains("_mean")) %>% + pivot_longer(!c(prop_NAs, model), + names_to = c("method"), + values_to = "auc_mean" + ), + auc_data %>% + select(c("prop_NAs", "model") | contains("_sd")) %>% + pivot_longer(!c(prop_NAs, model), + names_to = NULL, + values_to = "auc_sd" + ) %>% + ungroup() %>% + select(!c("prop_NAs", "model")) + ) %>% + mutate(method = method %>% + gsub( + pattern = "_mean", + replacement = "", fixed = TRUE + )) + # Preparing ARI data + ari_data <- averaged_data %>% select(c(prop_NAs, model) | contains("ari")) + + ari_data_long <- bind_cols( + ari_data %>% + select(c("prop_NAs", "model") | contains("_mean")) %>% + pivot_longer(!c(prop_NAs, model), + names_to = c("method"), + values_to = "ari_mean" + ), + ari_data %>% + select(c("prop_NAs", "model") | contains("_sd")) %>% + pivot_longer(!c(prop_NAs, model), + names_to = NULL, + values_to = "ari_sd" + ) %>% + ungroup() %>% + select(!c("prop_NAs", "model")) + ) %>% + separate(method, into = c("dim", "method", "useless"), sep = "_") %>% + select(-useless) + + return(list( + auc_data_long = auc_data_long, ari_data_long = ari_data_long, + auc_data = auc_data, ari_data = ari_data, + max_repetition = max_repetition, + sampling = sampling, + struct = struct + )) +} +``` +```{r plot-func, echo = FALSE} +plot_auc_ari_data <- function(df_list) { + auc_data_long <- df_list[["auc_data_long"]] + ari_data_long <- df_list[["ari_data_long"]] + max_repetition <- df_list[["max_repetition"]] + sampling <- df_list[["sampling"]] + struct <- df_list[["struct"]] + + auc_plot <- ggplot(auc_data_long) + + aes(x = prop_NAs, y = auc_mean) + + geom_line(aes(color = method)) + + geom_point(aes(color = method)) + + geom_ribbon(aes(ymin = auc_mean - auc_sd, ymax = auc_mean + auc_sd, fill = method), alpha = 0.2) + + ylim(c(min(auc_data_long[["auc_mean"]]), max(auc_data_long[["auc_mean"]]))) + + scale_x_continuous(breaks = scales::pretty_breaks(n = 10L)) + + ylab(TeX("\\bar{AUC}")) + + xlab("NA proportion") + + labs(fill = "Method", color = "Method") + + ggtitle(TeX(paste( + "$\\bar{AUC}\\pm s_{AUC}$", ", function of NA proportion. N=", max_repetition + ))) + + facet_grid(rows = vars(model), as.table = TRUE) + + theme_bw() + + + ari_plot <- ggplot(ari_data_long) + + aes(x = prop_NAs, y = ari_mean) + + geom_line(aes(color = method)) + + geom_point(aes(color = method)) + + geom_ribbon(aes(ymin = ari_mean - ari_sd, ymax = ari_mean + ari_sd, fill = method), alpha = 0.2) + + scale_x_continuous(breaks = scales::pretty_breaks(n = 10L)) + + ylab(TeX("$\\bar{ARI^d}$")) + + xlab("NA proportion") + + labs(fill = "Method", color = "Method") + + ggtitle(TeX(paste( + "$\\bar{ARI^d}\\pm s_{ARI^d}$", ", function of NA proportion. N=", max_repetition + ))) + + facet_grid(rows = vars(model), cols = vars(dim)) + + theme_bw() + + return((auc_plot | ari_plot) + + patchwork::plot_layout(guides = "collect") + + patchwork::plot_annotation( + title = "Metrics in function of NA proportions", + subtitle = paste0( + "With N = ", max_repetition, ", a ", struct, + " structure and a ", sampling, " sampling" + ), tag_levels = "A" + )) +} +``` + +```{r prep-fig-caps, echo = FALSE} +figcaps <- sapply(files_list, function(filename) { + data <- readRDS(file.path( + data_folder, + filename + )) + + max_repetition <- max(data[["repetition"]]) + sampling <- unique(data[["sampling"]]) + struct <- unique(data[["struct"]]) + + return(paste(struct, "structure with", sampling, "sampling and N =", max_repetition, "repetitions")) +}) +``` + +```{r plot, echo = FALSE, fig.cap=paste("Graph of metrics for NA robustness with a ", figcaps)} +for (filename in files_list) { + data <- readRDS(file.path( + data_folder, + filename + )) + df_list <- prepare_dataframes_auc_ari(data = data) + plot(plot_auc_ari_data(df_list = df_list)) +} +# filename <- files_list[1] +# data <- readRDS(file.path( +# data_folder, +# filename +# )) +# print(here()) +# df_list <- prepare_dataframes_auc_ari(data = data) +# knitr::knit_print(plot_auc_ari_data(df_list = df_list)) +``` \ No newline at end of file diff --git a/code/results/simulations/NA_robustness/NA_robustness_19-04-2024_15-18-55_uniform_nested_1-200.Rds b/code/results/simulations/NA_robustness/NA_robustness_19-04-2024_15-18-55_uniform_nested_1-200.Rds deleted file mode 100644 index 77041c7ad6d51f5a546aa167380c61f61cb9a561..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 12553 zcmY+J1yCHp7Os&1!8JH61PdMr65QRLpuwFD?z;pD?k*eL-CcqcWN~+QcNTbg_dcEa ztLv|x>2tcK&eV6R2N3Zd?mrjesUM0if!y?4eBKTf8bZ8Q0PUAx|4&Tep}DzuM7nfx zp9A^6V-*E_6U5&m|IECH5gve-0{KSRhtI+MmS*Da?(yFFe$C_jr3<*1^EDpZhbLF3FK3Nu#gW*8|0Y|KkY1M*@f2M^3{UI3E}V^iNZ?{U8+;|*AR%?t&f z86UnIe|rFPhD~0+t#mn~X1+NGkce5nOT#snj{TqJ4H%uxRYp;EY0h%X&Tt6k8-|U@ z-5D!8zuCqoQwxg)c(80Cbe2Mg0F!(E$T4TGd}A1@G9fQO6)1vhHl~kbD>)U>IM%|0 zD;AYbo*9uM8~o1bj;YFgj(@pf?TxZF|HT3$_|pZl+Ec9hvHh4bIvmA$Rh?VkN707g zJ6U>k=*3pzZYdM{{^Xk8$*~Nu`1#|;=As;pZWI+@zXvfE)4vOrif|#+PcKCy~-cqudbO zbP!(VKT~z>IrM}!^culygPmrtc)ap;2UFJ|}TD=bsI zTL~;C3@(u7| ztJVt^qwpIRCki4pNA0p%Dk;Ske1H5C@2h{iE09jqpy;Gg7^5DvFyyJQyhrAG4z;ic zV`YB7dLhI2zg9=4DT68X57*=PqTKGZipq z{QKx3pZCJuScH)HPBtINQrgaeTc*?6`}w7aJbSvxHj&*Et(XZ9*I^D2%UjtLSGMpV z925CwYmGe=)?ISns{xmmT>S;~M7`m(F4LIxV#p;j^9zVE){5j?5&T8;8e6ofBvEs$ z7uGUPwHS4ZVSQGKkJx1h#|4uqRT>r_#V$kwsdU5sG>}E?zmSpPOe3BufY zq;tUnj_B2iK*7#BqIq!n($KBS&lHdE1(h@BrrRc0rI5AnqpuB*mJ-#V@)VM+jzPFy zDC@K0bV8CEQTIe|qLmA&xhN_$a|-WZB5ww)GG&oYCHMhL@kM)YD$Mr|A!G z6drb^C!eq9$Xqm5@qaGV7MkC$y~P?~6-N&n>CaaATlN#7=eGDeI}(oZ`y#|*t8n(B z15LbB&etAx<7JJjWDQ(iW8klEn-6fdQ+Bn151tn##bSG+#Lw2j!W-uu86uR%kdh~jG=TYRN-bC|M^ zC?LMr#t8Y}^TMric&W0cPFQ;oYCw4f^e#1SP6ej5>p!Qv+Zw9APdrh7d?xgg_IvR} zp2$06j+9S;KZX71$=$N?ziY|J+b<73bh;Mc-uueqsJ2NhPYfc8Au2j=g4>zcopA_T z{pm6Ea{9%dqVU9L>3lDsnZ#dP(Ao7h740CeL}j|53yvR~prR}t z%!tXhFxXVOz#A$1&KSib;|C#5Z_-}xPOjvgQj6f((t@X`!%Ct&z~neP$#e_?T;O-M zD@`9~!DzZsJ~lw8@wyzh6svvI7rlYZ$h+;rJ2jFy!927!*xWHn4RXR)k1jd_2wG4z zhi`bm&vI=mI-M~y#aF17gXZS(ie+tuo6>wUdT z&aEmRAAcN#fC7LY7gW$WDgU8&dpjO~+bo?q0UcY(lvgL4;E$Yn2!qSkHYxfy5ss<; zAru;WW0u`uC;L%diXqB_Sk*LB^M2cv`@WBbngxJoz%S zF&qWNDO~TTr9GoAQyS-#0BRM$xhVW<*)z@Sd z*&TJb71{vX9LQuq?dqFOF-5+V};t^O3b&d zakuvUM{ab-i(!ja63)JnCmhaH7OfG(?^`HJpmHpLs%^zNC*0mfltrjfeut1;+L|e+ zpBxFVdLe*y9Wsidt^!sM2}xDzkUJB5J(v26fzxYA#TdAMYr_51f|DS|BKJn|QJY8C z9^(RpP4|2>VtZrk3$v5W2uWysiaW(7IxZ4X;M)5Nkh+?r#io5E*3MIjwBUZ4`l`M! z%V8*%fAX#tCGvF!PHPICra&o33DnP-`y7HX1-$Tg09)AJYSup3qex{ z*&=-<;{CI%Vxh#B2NbSkHQhQJhlu0)kH>AX{Cug(&xEedeXi1Vv@dw`eBVZeHSqS& zOrL_i*z%+Qr3UuZk*u)VdEPuIZ6TmL=UH}hhE+DKozwvnl#v8PHe#lh@Y9yRRby|SmX z??wD{XuGV9%NR${h*Yu0UjT;(XUI4MdVn+T`NdUr~0rexbd}? zt#kQzv)nJkaP3GL&)=m(*%487|Yu(ruMAzAwa6Ih<#U*9})U2&3-mAJ?goJ(m;2 z#mAb5RMB%GLG-H0+eVZM!bqE_N+S=WH-kcWAg}cC+@scWhi{9|afew)tc@h{*^CD( z;zao@HUXCL0ekhyk$tR0RaR@>LfdoRER}<~Q9?=mKJsODi|Pzk3QUKJ8n&0rR`?*E z33G9$k2djhta9mhd2&njRE*QEYG+gGU5a$iCbwlE*Cb7Id1H|_x#a04HRmx*VbIGF z+r(&mt|V|ddA9MV&YyOHFnqv1J`stEh6zJS&z;VXg}U@XEgnYK(uXb9lcy@4(~6T* z`$K((MOkHx40|037Gumef7Y)obbrU!nYsAtkwz%=>IM0Pf6^~^AbktY50X?isBK1> zj#`W+n_i~nA zGaaOJis#edFl8_i$EBA;J^7(b%aA!#^W}0YNV>nDfzG8z>0`BU`o<5L@1=^g&yH4G zM1c-$i}VXEh|6Yav@$O7iq+}~JDv%hgmRAV_0jx|6Kaeu>DdXXRKMLFb8=nJ*UI-G zs%&&=G)_++*k|P5V`h<)_TQN@1T#3}4MABCdu0b(FB4|>c8kfPW#kh@;HBeJcdD8_ z=M8DwY^QrLvpQXDpVbOjHCB^mk2HC6NAABj`11CT_{b-hO-IosGbA72RV^#jviaq$ z^(C*xU#(hSjj{bFIDAE4keae(sIjWV*}X&HCYASFyH24 zX2os5?Lj;_kofJLTnIwcQBZobpfqO?%pi_juHMRnBS%h-YlLfIx+ajQMqlA7EaJNg zL!3h(t{%82)CkkTo+x3SMn+bq*f?qDDC+U4$Co?aj%SbT{WaOL4eS73#V(=t#!E}) zP9ajYmNlFT&IJQj*q!QRamf)?3(n>*XHkb|W2fA&FM#*&eJW=a zp%$`i3U}l(J+X98KW5@3|Aa_N{VXq&{J6c4>W=sEbOl3vh5_*=U{B;j_v}}08dew= zoZWbX*l0n$Y)bqq5;qZrQrqiPg|e5Uh7_v-*xSZ~nh&%>gheAf+egzDgcZ5J-5WYc zFjOrWN8>3HEOQ8#Y)(KNF+;s#es3Pev(M@<=v^b(q;a{d0OU|~)$5i_rX}IPiw?dZ z>HPdA(m$TGcc@WECRo<)*0RQfU1+IIU^i@MyH+ zBm0L!X-$?iBe~>`Sfs_hOqL8HS>!lrhQn#?;#vBOBZ%GT%|p?!?!QjCOIbBRYvzeaI&B9ThxCBwkE#7Iq|Ii>`22vV`MGyR2X9`PUh zL!-Ms8@C%J5=%`F1=B_oPr9L9_a7^LL~Z)A$F6>orJrJ9eUIN>Nnn&hsKdUh_*)T? zd*e5F7%QA_>SOQvXY720@~d-SIs)tRbC?zr`VZ|Bhn=b}Jkef3Vc_Ij)`EvrX{+z!k`p2L@(J+N4$M?yEAp1DWS z+(7Vc7K$q`>XWq@aT96>C9C&x;@r3hMg3};qcMZrBfXgNK>Y6BjdOzbQQhn7XM1$- zX9(Rc>*20t#RBx{97SJ$^0Xl_vB}W^rVAFXT>7;}$(r7Z-jPWtv6^J0jGpPJksa^( z@)omQ;Ba{59)AVq`MOlu()nT=xbSyC!aA#C0t-{y4x2xtZbEEI9K}1PO(EF>=x}^} z4z6Zg`yw;^n{UeJj8;3k=wDulwd~jFDK^%}im#>>eavl*7KRRvl4%F#dWQ7;EsXSz z6=pq>3|}UP-^^W|vAqE{mqS9s&h*G6DsR|*RBsRy1H}c72vh;fp6Ik{Gy0R%|KD<$ z4F0>1?>-tt@#TM6meDt3Vt1BDc^?C77X9V|b@R#2tu~HDWz!fri*VqVZ%SQ+ZNp7^ z3chdLFA=rFuIBROkNrktcg2?x{lSpz@Zj}j)j^U{A3zfGYJ#T{|i5hmieU&-)eT87`2-@92`AL{O-A7TYbR;3t^p#jWv=9OgWvSN%etjC9?RjF(xV&+cE)N z(}wA=_Hd)STf}WP_`k!Es7e{w#n8@z^E!FBD8UZWcmYIpE1_o;Am8gX+RIf;J>pKx z$&Qh51YFCH7_0eC1kVx$yeFxM7Xz_xgXvx{R{!Yx#q>ER-fV55cr;-ZQ0e3p;D(<4Aw-;?i7A-kWl2_0LJ z#fZVeK=*;|@1iXbZ-Q|OxC88lRgrUIcyCp32}ddoZC1gdRLF(xp^y*&Q1jMQ2qVwEND{2x~CccLt3If4e^Yv9kA?G`wb|Z`- zhzW`nK{JB-fsg$*!zfUC_lH9VySRm0u#E4MXk4dKojNN7iop2=c#lS+2xRZrVx8Pf zU2{NDu);>!$uRA$nO%qTa$^|WTCovja@A_o^@TT2m>R&38|W4QHrwp$*JNF%4zb!l zlp?>sw++0%I@1LUeyMvPex)#Eac5aVF1Gmj)agP@N!m*Gr3=z~k`Vs;dKf--%r7sx zEDa&1%iN`-#M<}qXrHU7JOTdoeMK+k5xu|B1VGol4|$?>#qIB~foMdlyQ=Rc@x&Mr zA>oI6WLHi1DeZlr{~bl*t&2WbRNWmTLl2~DzHS#d`2#3Gh14$-Qf3BfbvRk9_z*|@ zn!+5i8U;Vo_=$f(7VCPjQ$VNf-_f0a%6@~yJ$&5a)(ZGL)&;Wl3BSIRTE2?SqGuNg z(t;pISUIyw4~*M-CcQ4C=8;tB?LpsOP%BwzlAMP)mbW9~UBIFLUq}6}nGVX>8*H(J zCIh&&{pl~}$Q3w;`fc@Z(A(estrts0CuHact+O_0$13ngEkD#}#x49yWlq=PA~u=8 zuuEG2K?ZvOXC2L+|AUd6YfS?3+Rx&cr!KjOn~JRG$L`Nq*$h@`5T$wKdS9A&-gV-M z-e1QkdsqXo=9sZ;vors_{EiZy8YC|24ZGqw2GLH-5>FzEo;~0~Bw|%#_)W0V`}yUc zvtL-VevoA82>*06>?}u^`FxRC<_l(7cYSrI-1E@Es5)tDAc>Wr#8NkALo{`<7;8+z zQuQ*^T#T9NFcqoVI?h5M2DR3DygO%x8%B4X(W@i67TJE6m$nUmdMQj>OKjKle7j(@ zCQ^N?;1FH2=p^b+mp*phO6n;>(J>aU*bj0Xo-fOWN{hxoP_=58LMp(4h5 zA?X=Eqee{Bj6p?NYI8@SV7QCBq{41#y) z>{p*`K`C|W`RvMJ_Zdre&@luv5duk+BGlm1wJR@hrqE_m1rg4;lDbmM$J6 z!i>Rfnz{6>?zum)8OUx~f8t`Q>o7*twUNGCHTC`VAgjWuFRd4%5|9>Q9 zFp6#%fAw!**$b|uANb7f%JIL&ZN>ZV;eRS45!AEG^2$+g<+hTrBh&pC3Plb*ji|NQ znW4s6@RWrS+CDM6)aPu`z+aO})cQ|LJ2Gb#W&{7tvAC+KW_H@o14B#YjhOBcVEGvtZ0A+QAAYk z*;_jxjs#H%$m$IRQE2Qco9=opeN+F9D&Xi~+ES|r(ftp4SLj8YffQf`+@Qpqw1NpZN{o2%OGZn)B5dE#J zH<}I*wS!8Y2N)e{gzlDtgwvoOaCA(jStuhg^HJeGmLtDZDg% zK83ack!+=eIRTjLcz)NmrHd7`40;obJa7bTu?K#Eb?vYEH`%?Ini`h}lHNBdXyoIU z+$_66ZovGx-Bmyx(D5m+=r4oc8%8ULOsu~nFL;}fE5|cgX7kY?4N|T;=R2!0H0c~DqbOXn4czF+s?X2F+OoUCTH>vy+f16?5w61cxPk| z8JYIYlHQcz=bE`YQx&lhrwFDA&Vo$c|E5g&L>}C3KR7vND38@U4T{U!W0}O?)ZrWY zn>riTqm&BOCn;?f24Y|ASjKp-js23>c_M`olk*`mg(DMHlJKI3vtf}xUtlkbL}kV5 zlwU`Yy$Vq>MIyG%Gms!QUFSgBZVA;NwB~!v^bvViUhC=I35Yp0XLS)M`LfD{kPz+vEC_7I0hKsJ_EBfWBAlTGW+7>!4Ds_ zQPwvVT?5*~4%x8f@70@8UzC5WtJ|etze5ud+rM@86Sv~|qc${x8BU229v(7pX4?Ik z-Z%&l^xVba`EIrg#YTp1Kos-bpkcT(m`e8f#Vy5$OzOdbkLXK}$giv3@Ym=@NN3w0 zg)!{)xwDbik?Rz@H|6buc!I>il`~>5EfS-ZDZL!hm!CZ<2Oe6;cH1k@K$Ol|Y~Z=m zt`baT$k>Oz zQv088ayeSWGfc!y?6aMPafjnKs??e_GZ=7B^-&!3JE0an*sQ;LQSWbxjQC?1FF9m^ zQuTssZW#>yB@$&!XXU)nnq8uqo}tfkGHT9@CPNJpH?sH;S8%CHR9cV3``0qlffE>f zoIUZjcEYy>GDps8^+P|WX?U6jfG6a;A6_Q!(GhYihXoctMq|u;Iz@99TtRtQEcMu6^;)%9ZnBX~zXiwFc3*Ilbb;SXOW#n6Oj@ zmfFxXFT^wk1yp1b@{pC+Ch_+LEl`VKyzqTh#-#y!?QYqgB z4+TVZitO*W<9JrC*%&&LzB~KwO6P~7iFsk}qx16aRqI6!WljirZQ}MmdLJ6}S7oK8 zi6`3sZ+Eo-lyqEl_soS9WWntJTvOJ628Z-%t_9Tm68w`2C-F2vekU-ybXBH>%Kr#6 z?Jp;C1P7iu50)1 zUST8Th6Ed}tI3ZLBdeqcp#kxc%-h}vDFC}>4&Cr5+tx4(>VsDWe;?Cu@8Bv751bwy zL~C=&$=$ceOW6;|R^+#T84On0aGwi0PeR#K^PDBdBmta+jcrfTpR6htsN}5kswa-$pqs*>*w<Q(i;C)mlmOer*Nf zvZ2h`*!e_Hpt6kAepr@eaBiirISAoA_^Zxt_NmG3WzOsz^HXcdcUO-|6SSt%uh*KZ zk&teYnP7cgRcvmd5sCpPARPE=GERgU3jsLyi^*bJ9!=$}-YUewY#n|cxJ=o4err`dlj{5bcHnD$aIXu@Mn%b8RD!%aVW8r-Lpt>FeO5u>yr;dJR9l`Fu6a zd8&1Hy4WZy`jxV0j|-x(B+{NI6puO@3O>}M$gy3|5%YOCq`C+T^lZq>eQYhdJGj4v zInnPU#b;NgyfiQoIelYpx00s#iEay#{>yF~^J8fXVHadl9vI`6<2%0{h^PS-{d1&( z9$(Ka`4AiPsKNcgYvL;CU})FuN;z=xvNvX9bbp;K@T_OvU2$OR!uHnoN>M~)466qBu_nJ zBabs01m2lEC}$uS&G%drul`}iUQm6u+I$UOQb}mHmG^x~t2edH>Q&LJcPmvU$WG9X z*7rW=zi5?aAJUs8GD(os@>)tTesu1GU(2SpG%|sSJ(s$>O?x=k_W5=F3S-qk266Mtj=T zcT4O>S@OzL@aHm_??;xB&-Xolh>Cd|_RphF_FXvcj53?(iEF^7&G17#&yt!*j~~hm zp4zIWrWXaxyVpYvN?$adlSSRyE4P_;iP4};tA@`M&C=b=oXwz+eXUPjRv;ORaRtw` z9O8$|Y_#Bn9{{CI*)x~Y^MV^7@ipWpWIFB7Fhw_*7u6QeZYwD)0`%W zt@9f4d8|{tw0B>Ww6OUEZdrcbYC6n)*?^|MCbVJ~j#G<|XC?Z=uE$jT$FIGbEPE+F zRWWh2+`JQ*Ph36K<;)jp$7q>&Z6Bp;K7VZA86c`F+{eW$S~BV_ud-SXAzl4^YN)-R zEw98cR6LHMT4I;9@}lq0?m~A~kjSWMPxfl+h$8hV^K}%Y+aZqsmX{?{Dtw~k`dHKVzilxE} zTg{^Hf6pHni0Ys=B^t-nxy3@Tg(xxqX0H}5D(&uNYrg~2lq^;4%lfmQAz}%)|60I+ zCxUAad+s9_lCx*(7Z>>_(cLgy5d#3B_&0UeOE-@b-}3aN$GPtxRt}n(nfJejp3yuz zpN-D42R1#A@!Hty1e5l%?{A)o*@6vsX&q|HP0%hp(pv#XzaBT4wJLlj=Tu94nTyVg z_)iRc6~g?U@C321--K_~LajxvUB_2Xu?57UaxW!1Q214~t{nUpE{12kf~vlB>w0j_ zf&K24DF`v&Ii|MfEwy^7I@Qge{z2#I$;d80+n#?A9lIJkpMolAWG)$tU|${I=JRU0 zW!ik>MP1&f9siwbnS?atHb=hiF6}?ctbg6$iKU5#MRCYYOD11L&)V+gC9L_j=THNh z!pFIRhGrT6Iz>N50~-5W`Y);EDVL-}J=_r1u3n=_G%8JBDb9WG+Zhst{Qi11aN`-h zS*o=hPP#UmdC|<79>L>+Zmx&c^aHU~DudH$+8_1fp@>}z_DsNi)(x`NGPHMq4Mymq zkWkVi1{!EcMfUBv-f_tJwX` z6c4!7T1nn~c>vSd=_P+fr08Q_1dU6~kIlyz?56oUOO0mK`@o*-_adAfTa&S0?YqWD zJl+TpwC40lee>NB_*cwBj`uhHk-=<7e*l>6y;jMu$pYB&h%tu)K}Vq%H22&$*YAOB zOHO%zX*Zsjav==qVrq*uK&Ufjj_dc0Y6jwQC@-|cF69T|*{7}wOrw3;o|rAhmWEz$ zlh^Bzrd3&ddx)DnD4QGg!O+V0$(|r|RVF^!!yxn&yH{~EzSVTUkgf8Hr2tj+qsQ*n zcy7T**e}pfM)v0Wfb?rZ_@v?7YuU7KQqMYJ-BsM&OpOM=HaC)Wxbxs+Pa?xzP(!ju zshT!sDJ48t8@5|w(qg0@;+oaIo>i_T>I78;>ZnH)8vHo};49KVQm%diJXoNGH?;pa zQmk}MIs?K_OlEUiF4vOXWDty@I>v7jnfx>G$FHya1$cRBSaRTQ!a+5)+Hz}V`qQV4 z4E%WSrMmW|*e=_WM^Ir!gloPL?`*eey)9rmu&1JMyzSxO14dK3<}JUU9ZsZSQYW(u z>-$&h>MHc!IURNyK+>Xue^Y(X;AXJcZV8>OLPD0<;~pMw8++a6Ev=*ROkwpfsMY(8 zZ}CFRxVXNX_u-E-A3M&Bun@=RLWGGN5}PaIHZi)979W`pI-V5UeL!L-0wXHqWrx@AGXGZb^J6=D$o1fSyPM?&icX^U5<4LYzr9e@W;NwMfU^{8D$~z6OJ3U z#Hgs@R_$li9q_;Uz~iH4&qUYy8`0=>8-`eq4PbROR)Iu}1a^$u8wZZPU9DJo^I~j- zN%V|6A?9sr+bm=rdgOnC1?}V7u6n0jI#yOuYldV}fNHMgaO_<@1#Ms7;h3-3RB7#F z9n_FGN(c#yrUA#RuM^J1H*rD~&2lYU3h4owAUR7t2Y)lYB4) z&wi9WcR=Ln@|@w{YMm^&B=|b4hc=^NB1U;S=hu+*a_}tmVkupc+5*U0s}c5@7Y75R z-DThMkxF?b?J8##Isaer>#PT2ik_-a)~ z(jq{CQctDvtQLz&#d*KIvNS9a>|_Js;fXQIyqwgwe(#UqQaL|49DcE&Rw+BauK%#T zc(#p_x&;+leW_yFzQjDsyC2$EKii}%wDcF*cL&`2zN4f(r6j+34?fuv<{V=D8Qb;E zUj)$|+raOLH#4S*@0w(vTA0^vZ|iYU2#Wv6tZ`x0be1EQR^s{EjOd*a_QXJlcRSeO z8!`B~IpqAr4UZi=87e#F*q~Rt zTGOYV=x?tvVON79m(2+Gsrtb*m#*ULX-7dT3lu*kMMx9(aNNrco0DbPm3y-u+hP`d z05Nhd_opv#F>g+QWiiHtgQ8_OLBbxtT>DX{Ol%2gA|CAW9tQ#X$cu+cE+K&P@D94S zN^|j}gBt}1-shXkz%ci9@S^?ZJy`|~LM>nk?JU3$04 zDAg=QU$g5fu5#r5sW;ij^`oX1$yar%0Xxl61Rv9tlOp&2W0#!(w#`GvbJnRJu-sA= zoL5|k%U2S1D-l_P$ZK=D5|wuzEtFy#LdUcG2B6x)J5_|PTQi-O?d{wZn|B;7 zAg-YUT=#Ybh5>=ttmZgbloL04*c4OfMK3DzARdU+Jw{=J^s1isXH3`IaIa)5i{!8B zwjV77KSR#*ymO7)@K8n2FV7#4fV{H(ok7GQi_;0H9r&1q_`6e1hwqH)dgpLn6JOd9 zC!uI|J?t{;@608Xtc=ij#_XR;DH{Xeaj{Obo~jzh%}Wc0hf z^$C9NB9dC*PVOjoq9oSGdC=;h9q`(DuLbN8x#?oWyYs7AkrQ}G1gTZdRKKgHt-mm#(*E{^qE4s&j%CY z79=yY{6E31af+_drB0Z&`d{QZ;j>3F zz%JxTI6!}pbUJAr2}&g@kV$wC;@_sq8FW`hHUVC^mckD?y%iS!%eG9}-S z)~10Ss9Pb#nAEZJZcS!K3|2H;Oc-$l;FF|{hFS0?*Gb}bIA!k8rdQkJ36Mkc21`Go zC(+RkFucP)X2;Nqjfjf*j*f_zh`)$>VjiuV8Ypfp?|Gq$ZOmB`&d|GUCis=6L~wUP zQrQ&7#la|uf{!<%#$aKsH7h5k*T^J5&qWkxbdz}5?a(&Vr%Em|2Yo}C4!_VWd4hxcA9;={`2YX_