diff --git a/DESCRIPTION b/DESCRIPTION index 9d6cc98..86f71cc 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,8 +1,8 @@ Package: discord Type: Package Title: Functions for Discordant Kinship Modeling -Version: 1.0.1 -Date: 2021-03-19 +Version: 2.0.0 +Date: 2021-05-15 Authors@R: c(person("S. Mason", "Garrison", email = "garrissm@wfu.edu", role = c("aut", "cre")), person("Jonathan", "Trattner", email = "code@jdtrat.com", diff --git a/R/func_discord_data.R b/R/func_discord_data.R index 44336a9..2ca1951 100644 --- a/R/func_discord_data.R +++ b/R/func_discord_data.R @@ -8,11 +8,9 @@ #' @param id A unique kinship pair identifier. #' @param sex A character string for the sex column name. #' @param race A character string for the race column name. -#' @param pair_identifiers A character vector of length two that contains the -#' variable identifier for each kinship p -#' @param demographics Indicator variable for if the data has the sex and race -#' demographics. If both are present (default, and recommended), value should -#' be "both". Other options include "sex", "race", or "none". +#' @param pair_identifiers A character vector of length two that contains the variable identifier for each kinship p +#' @param demographics Indicator variable for if the data has the sex and race demographics. If both are present (default, and recommended), value should be "both". Other options include "sex", "race", or "none". +#' @param legacy Logical Logical: FALSE (by default) when true uses legacy code version #' #' @return A data frame that #' @@ -28,7 +26,17 @@ #' race = NULL, #' demographics = "none") #' -discord_data <- function(data, outcome, predictors, id = "extended_id", sex = "sex", race = "race", pair_identifiers, demographics = "both") { +discord_data <- function(data, + outcome, + predictors, + id = "extended_id", + sex = "sex", + race = "race", + pair_identifiers= c("_s1", "_s2"), + demographics = "both", + legacy=FALSE, + ...) { +if(!legacy){ # non-legacy version #combine outcome and predictors for manipulating the data variables <- c(outcome, predictors) @@ -53,12 +61,129 @@ discord_data <- function(data, outcome, predictors, id = "extended_id", sex = "s if (demographics == "none") { output <- out %>% purrr::reduce(dplyr::left_join, by = c("id")) } else if (demographics == "race") { - output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(race, "_s1"), paste0(race, "_s2"))) + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(race, pair_identifiers[1]), + paste0(race, pair_identifiers[2]))) } else if (demographics == "sex") { - output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(sex, "_s1"), paste0(sex, "_s2"))) + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(sex, pair_identifiers[1]), + paste0(sex, pair_identifiers[2]))) } else if (demographics == "both") { - output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(sex, "_s1"), paste0(sex, "_s2"), paste0(race, "_s1"), paste0(race, "_s2"))) + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", + paste0(sex, pair_identifiers[1]), + paste0(sex, pair_identifiers[2]), + paste0(race, pair_identifiers[1]), + paste0(race, pair_identifiers[2]))) } + }else{ + arguments <- as.list(match.call()) + y <- ysort <- NULL + + IVlist <- list() + outcome1=subset(df, select=paste0(arguments$outcome,sep,"1"))[,1] + outcome2=subset(df, select=paste0(arguments$outcome,sep,"2"))[,1] + + #create id if not supplied + if(is.null(id)) + { + id<-rep(1:length(outcome1[,1]))} + #If no predictors selected, grab all variables not listed as outcome, and contain sep 1 or sep 2 + if(is.null(predictors)){ + predictors<-setdiff(unique(gsub(paste0(sep,"1|",sep,"2"),"",grep(paste0(sep,"1|",sep,"2"),names(df),value = TRUE))),paste0(arguments$outcome)) + #unpaired.predictors=setdiff(grep(paste0(sep,"1|",sep,"2"),names(df),value = TRUE,invert=TRUE),paste0(arguments$id)) + } + + + if(!doubleentered){ + outcome2x<-outcome2 + outcome2<-c(outcome2[,1],outcome1[,1]) + outcome1<-c(outcome1[,1],outcome2x[,1]) + + if(scale&is.numeric(outcome1)){ + outcome1<-scale(outcome1) + outcome2<-scale(outcome2) + } + DV<-data.frame(outcome1,outcome2) + DV$outcome_diff<- DV$outcome1-DV$outcome2 + DV$outcome_mean<-(DV$outcome1+DV$outcome2)/2 + + remove(outcome1);remove(outcome2x);remove(outcome2) + + for(i in 1:length(predictors)){ + + predictor1x= predictor1=subset(df, select=paste0(predictors[i],sep,"1"))[,1] + predictor2=subset(df, select=paste0(predictors[i],sep,"2"))[,1] + predictor1<-c(predictor1[,1],predictor2[,1]) + predictor2<-c(predictor2[,1],predictor1x[,1]) + if(scale&is.numeric(predictor1)){ + predictor1<-scale(predictor1) + predictor2<-scale(predictor2) + } + remove(predictor1x) + IVi<-data.frame(predictor1,predictor2) + IVi$predictor_diff<-IVi$predictor1-IVi$predictor2 + IVi$predictor_mean<-(IVi$predictor1+IVi$predictor2)/2 + names(IVi)<-c(paste0(predictors[i],"_1"),paste0(predictors[i],"_2"),paste0(predictors[i],"_diff"),paste0(predictors[i],"_mean")) + IVlist[[i]] <- IVi + + names(IVlist)[i]<-paste0("") + } + }else{ + + if(scale&is.numeric(outcome1)) + + {outcome1<-scale(outcome1) + outcome2<-scale(outcome2) + } + DV<-data.frame(outcome1,outcome2) + + DV$outcome_diff<-DV$outcome1-DV$outcome2 + DV$outcome_mean<-(DV$outcome1+DV$outcome2)/2 + + remove(outcome1);remove(outcome2) + for(i in 1:length(predictors)){ + predictor1=subset(df, select=paste0(predictors[i],sep,"1"))[,1] + predictor2=subset(df, select=paste0(predictors[i],sep,"2"))[,1] + if(scale&is.numeric(predictor1)) + {predictor1<-scale(predictor1) + predictor2<-scale(predictor2) + } + IVi<-data.frame(predictor1,predictor2) + IVi$predictor_diff<-IVi$predictor1-IVi$predictor2 + IVi$predictor_mean<-(IVi$predictor1+IVi$predictor2)/2 + names(IVi)<-c(paste0(predictors[i],"_1"),paste0(predictors[i],"_2"),paste0(predictors[i],"_diff"),paste0(predictors[i],"_mean")) + IVlist[[i]] <- IVi + names(IVlist)[i]<-paste0("") + } + } + + + DV$id<-id + DV$ysort<-0 + DV$ysort[DV$outcome_diff>0&!is.na(DV$outcome_diff)]<-1 + + # randomly select for sorting on identical outcomes + + if(length(unique(DV$id[DV$outcome_diff==0]))>0){ + select<-sample(c(0,1), replace=TRUE, size=length(unique(DV$id[DV$outcome_diff==0&!is.na(DV$outcome_diff)]))) + DV$ysort[DV$outcome_diff==0&!is.na(DV$outcome_diff)]<-c(select,abs(select-1)) + + } + DV$id<-NULL + names(DV)<-c(paste0(arguments$outcome,"_1"),paste0(arguments$outcome,"_2"),paste0(arguments$outcome,"_diff"),paste0(arguments$outcome,"_mean"),"ysort") + + merged.data.frame =data.frame(id,DV,IVlist) + + id<-ysort<-NULL #appeases R CMD check + + merged.data.frame<-subset(merged.data.frame,ysort==1) + merged.data.frame$ysort<-NULL + merged.data.frame <- merged.data.frame[order(merged.data.frame$id),] + if(!full) + {varskeep<-c("id",paste0(arguments$outcome,"_diff"),paste0(arguments$outcome,"_mean"),paste0(predictors,"_diff"),paste0(predictors,"_mean")) + + merged.data.frame<-merged.data.frame[varskeep] + } + output<-merged.data.frame + } return(output) diff --git a/R/func_discord_regression.R b/R/func_discord_regression.R index 1ad1207..d7303c8 100644 --- a/R/func_discord_regression.R +++ b/R/func_discord_regression.R @@ -10,6 +10,7 @@ #' @param race A character string for the race column name. #' @param pair_identifiers A character vector of length two that contains the variable identifier for each kinship pair. #' @param abridged_output Logical: FALSE (by default) and the fit model will be summarized with the \link[broom]{tidy} function. FALSE and the full model object will be returned. +#' @param legacy Logical Logical: FALSE (by default) when true uses legacy code version #' #' @return Either a tidy data frame containing the model metrics or the full model object will be returned. See examples. #' @@ -35,13 +36,17 @@ #' abridged_output = FALSE) #' discord_regression <- function(data, - outcome, - predictors, - id = "extended_id", - sex = "sex", - race = "race", - pair_identifiers = c("_s1", "_s2"), - abridged_output = FALSE) { + outcome, + predictors, + id = "extended_id", + sex = "sex", + race = "race", + pair_identifiers = c("_s1", "_s2"), + abridged_output = FALSE, + legacy=FALSE, + ...) { + +if(!legacy){ # non-legacy version check_discord_errors(data = data, id = id, sex = sex, race = race, pair_identifiers = pair_identifiers) @@ -90,6 +95,26 @@ discord_regression <- function(data, model <- model %>% broom::tidy() } +}else{ + if(!discord_data){ + data<- discord_data(outcome=outcome,doubleentered=doubleentered, + sep=sep, + scale=scale, + data=data, + id=id, + full=FALSE, + legacy=TRUE) + } + arguments <- as.list(match.call()) + if(is.null(predictors)){ + predictors<-setdiff(unique(gsub("_1|_2|_diff|_mean|id","",names(data))),paste0(arguments$outcome)) + } + if(is.null(additional_formula)){ + additional_formula="" + } + model<-lm(as.formula(paste0(paste0(arguments$outcome,"_diff"," ~ "),paste0(predictors,'_diff+',collapse=""),paste0(predictors,'_mean+',collapse=""),arguments$outcome,"_mean",paste0(additional_formula))),data=data) + +} return(model) diff --git a/hidden/func_discord_data_alt.R b/hidden/func_discord_data_alt.R new file mode 100644 index 0000000..44336a9 --- /dev/null +++ b/hidden/func_discord_data_alt.R @@ -0,0 +1,65 @@ +#' Restructure Data to Determine Kinship Differences +#' +#' @param data A data frame. +#' @param outcome A character string containing the outcome variable of +#' interest. +#' @param predictors A character vector containing the column names for +#' predicting the outcome. +#' @param id A unique kinship pair identifier. +#' @param sex A character string for the sex column name. +#' @param race A character string for the race column name. +#' @param pair_identifiers A character vector of length two that contains the +#' variable identifier for each kinship p +#' @param demographics Indicator variable for if the data has the sex and race +#' demographics. If both are present (default, and recommended), value should +#' be "both". Other options include "sex", "race", or "none". +#' +#' @return A data frame that +#' +#' @export +#' +#' @examples +#' +#' discord_data(data = sample_data, +#' outcome = "height", +#' predictors = "weight", +#' pair_identifiers = c("_s1", "_s2"), +#' sex = NULL, +#' race = NULL, +#' demographics = "none") +#' +discord_data <- function(data, outcome, predictors, id = "extended_id", sex = "sex", race = "race", pair_identifiers, demographics = "both") { + #combine outcome and predictors for manipulating the data + variables <- c(outcome, predictors) + + #order the data on outcome + orderedOnOutcome <- purrr::map_df(.x = 1:base::nrow(data), ~check_sibling_order(data = data, + outcome = outcome, + pair_identifiers = pair_identifiers, + row = .x)) + + out <- NULL + for (i in 1:base::length(variables)) { + out[[i]] <- purrr::map_df(.x = 1:base::nrow(orderedOnOutcome), ~make_mean_diffs(data = orderedOnOutcome, + id = id, + sex = sex, + race = race, + pair_identifiers = pair_identifiers, + demographics = demographics, + variables[i], row = .x)) + } + + + if (demographics == "none") { + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id")) + } else if (demographics == "race") { + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(race, "_s1"), paste0(race, "_s2"))) + } else if (demographics == "sex") { + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(sex, "_s1"), paste0(sex, "_s2"))) + } else if (demographics == "both") { + output <- out %>% purrr::reduce(dplyr::left_join, by = c("id", paste0(sex, "_s1"), paste0(sex, "_s2"), paste0(race, "_s1"), paste0(race, "_s2"))) + } + + return(output) + +} diff --git a/hidden/func_discord_regression_alt.R b/hidden/func_discord_regression_alt.R new file mode 100644 index 0000000..4ce0556 --- /dev/null +++ b/hidden/func_discord_regression_alt.R @@ -0,0 +1,89 @@ +#' Perform a Linear Regression within the Discordant Kinship Framework +#' +#' @param data A data frame. +#' @param outcome A character string containing the outcome variable of +#' interest. +#' @param predictors A character vector containing the column names for +#' predicting the outcome. +#' @param id A unique kinship pair identifier. +#' @param sex A character string for the sex column name. +#' @param race A character string for the race column name. +#' @param pair_identifiers A character vector of length two that contains the variable identifier for each kinship pair. +#' @param abridged_output Logical: TRUE (by default) and the fit model will be summarized with the \link[broom]{tidy} function. FALSE and the full model object will be returned. +#' +#' @return Either a tidy data frame containing the model metrics or the full model object will be returned. See examples. +#' +#' @export +#' +#' @examples +#' +#' # Return an abridged model output using the \link[broom]{package}. +#' discord_regression(data = sample_data, +#' outcome = "height", +#' predictors = "weight", +#' pair_identifiers = c("_s1", "_s2"), +#' sex = NULL, +#' race = NULL) +#' +#' # Return the full model output. +#' discord_regression(data = sample_data, +#' outcome = "height", +#' predictors = "weight", +#' pair_identifiers = c("_s1", "_s2"), +#' sex = NULL, +#' race = NULL, +#' abridged_output = FALSE) +#' +discord_regression <- function(data, outcome, predictors, id = "extended_id", sex = "sex", race = "race", pair_identifiers = c("_s1", "_s2"), abridged_output = TRUE) { + + check_discord_errors(data = data, id = id, sex = sex, race = race, pair_identifiers = pair_identifiers) + + if (is.null(sex) & is.null(race)) { + demographics <- "none" + } else if (is.null(sex) & !is.null(race)) { + demographics <- "race" + } else if (!is.null(sex) & is.null(race)) { + demographics <- "sex" + } else if (!is.null(sex) & !is.null(race)) { + demographics <- "both" + } + + preppedData <- discord_data(data = data, + outcome = outcome, + predictors = predictors, + id = id, + sex = sex, + race = race, + pair_identifiers = pair_identifiers, + demographics = demographics) + + # Run the discord regression + realOutcome <- base::paste0(outcome, "_diff") + predOutcome <- base::paste0(outcome, "_mean") + pred_diff <- base::paste0(predictors, "_diff", collapse = " + ") + pred_mean <- base::paste0(predictors, "_mean", collapse = " + ") + + + if (demographics == "none") { + preds <- base::paste0(predOutcome, " + ", pred_diff, " + ", pred_mean) + } else if (demographics == "race") { + demographic_controls <- base::paste0(race, "_s1") + preds <- base::paste0(predOutcome, " + ", pred_diff, " + ", pred_mean, " + ", demographic_controls) + } else if (demographics == "sex") { + demographic_controls <- base::paste0(sex, "_s1 + ", sex, "_s2") + preds <- base::paste0(predOutcome, " + ", pred_diff, " + ", pred_mean, " + ", demographic_controls) + } else if (demographics == "both") { + demographic_controls <- base::paste0(sex, "_s1 + ", race, "_s1 + ", sex, "_s2") + preds <- base::paste0(predOutcome, " + ", pred_diff, " + ", pred_mean, " + ", demographic_controls) + } + + model <- stats::lm(stats::as.formula(paste(realOutcome, preds, sep = " ~ ")), data = preppedData) + + if (abridged_output) { + model <- model %>% + broom::tidy() + } + + return(model) + +} diff --git a/hidden/sysdata.rda b/hidden/sysdata.rda new file mode 100644 index 0000000..a68e1f2 Binary files /dev/null and b/hidden/sysdata.rda differ diff --git a/man/discord_data.Rd b/man/discord_data.Rd index 5329e11..320c1a5 100644 --- a/man/discord_data.Rd +++ b/man/discord_data.Rd @@ -11,8 +11,10 @@ discord_data( id = "extended_id", sex = "sex", race = "race", - pair_identifiers, - demographics = "both" + pair_identifiers = c("_s1", "_s2"), + demographics = "both", + legacy = FALSE, + ... ) } \arguments{ @@ -30,12 +32,11 @@ predicting the outcome.} \item{race}{A character string for the race column name.} -\item{pair_identifiers}{A character vector of length two that contains the -variable identifier for each kinship p} +\item{pair_identifiers}{A character vector of length two that contains the variable identifier for each kinship p} -\item{demographics}{Indicator variable for if the data has the sex and race -demographics. If both are present (default, and recommended), value should -be "both". Other options include "sex", "race", or "none".} +\item{demographics}{Indicator variable for if the data has the sex and race demographics. If both are present (default, and recommended), value should be "both". Other options include "sex", "race", or "none".} + +\item{legacy}{Logical Logical: FALSE (by default) when true uses legacy code version} } \value{ A data frame that diff --git a/man/discord_regression.Rd b/man/discord_regression.Rd index b61d860..2debff2 100644 --- a/man/discord_regression.Rd +++ b/man/discord_regression.Rd @@ -11,7 +11,10 @@ discord_regression( id = "extended_id", sex = "sex", race = "race", - pair_identifiers = c("_s1", "_s2") + pair_identifiers = c("_s1", "_s2"), + abridged_output = FALSE, + legacy = FALSE, + ... ) } \arguments{ @@ -30,16 +33,20 @@ predicting the outcome.} \item{race}{A character string for the race column name.} \item{pair_identifiers}{A character vector of length two that contains the variable identifier for each kinship pair.} + +\item{abridged_output}{Logical: FALSE (by default) and the fit model will be summarized with the \link[broom]{tidy} function. FALSE and the full model object will be returned.} + +\item{legacy}{Logical Logical: FALSE (by default) when true uses legacy code version} } \value{ -A tidy dataframe containing the model metrics via the - \link[broom]{tidy} function. +Either a tidy data frame containing the model metrics or the full model object will be returned. See examples. } \description{ Perform a Linear Regression within the Discordant Kinship Framework } \examples{ +# Return an abridged model output using the \link[broom]{package}. discord_regression(data = sample_data, outcome = "height", predictors = "weight", @@ -47,4 +54,13 @@ pair_identifiers = c("_s1", "_s2"), sex = NULL, race = NULL) +# Return the full model output. +discord_regression(data = sample_data, +outcome = "height", +predictors = "weight", +pair_identifiers = c("_s1", "_s2"), +sex = NULL, +race = NULL, +abridged_output = FALSE) + } diff --git a/testing.Rmd b/testing.Rmd deleted file mode 100644 index 38582f0..0000000 --- a/testing.Rmd +++ /dev/null @@ -1,153 +0,0 @@ ---- -title: "R Notebook" -output: html_notebook ---- - -```{r setup} -library(tidyverse) -library(discord) -library(tictoc) -``` - -```{r generate simulation data} - -# this function allows you to generate any paired data for any level of relatedness. -# r_all is the relatedness coefficient. -# Default is MZ twin vs DZ twin. -# -# in this, y1_1 is one variable for sibling one, y1_2 is one variable for sibling 2. -# y2_1 is one variable for sibling one, and y2_2 is one variable for sibling 2. - -SIMDATA <- discord:::kinsim_multi(r_all = 1, - npg_all = 1200) %>% - tibble() %>% - relocate(c(id, r), .before = A1_1) - -``` - -```{r test on Mason function} -set.seed(18) -MasonRegression_MZTwins <- discord:::discordDataUpdating(SIMDATA, outcome = "y1", predictors = "y2", id = "id", - sex = NULL, race = NULL, pair_identifiers = c("_1", "_2"), - demographics = "none") %>% -discord:::discord_regression(predictors = "y2", outcome = "y1") %>% broom::tidy() - -``` - -```{r examine results of Mason function} - -# if this works, we should see not see a significant difference score since the covariance of the ACE are -# 0 by default. If we wanted to specify cov_a = 1, cov_c = 1, and cov_e = 1, then we would find a significant -# difference score. The two variables would be very highly correlated. - -# cov_a is the covariance for a1 and a2 (variable 1 and variable 2's added variance) -- -# how much genetic component overlaps -# this translates into the genetic aspect of the correlation -# -# A and C are the familial covariance. We cook out the variance associated with A and C, and the only thing that should -# signal a significant difference score is the covariance of E. - -# if we want to be super confident, generate 200 datasets. Write a function to flag whether it's significant or not, and then -# count proportion. -MasonRegression_MZTwins - -``` - -```{r try my function} - -set.seed(18) -JTRegression_MZTwins <- discordRegressionUpdating(data = SIMDATA, - outcome = "y1", - predictors = "y2", - sex = NULL, - race = NULL, - pair_identifiers = c("_1", "_2"), - id = "id") - -``` - -```{r compare our functions} -waldo::compare(MasonRegression_MZTwins, JTRegression_MZTwins) -``` - -```{r define significants function} - -isSig <- function(df) { - - model <- discord_regression(data = df, - outcome = "y1", - predictors = "y2", - sex = NULL, - race = NULL, - pair_identifiers = c("_1", "_2"), - id = "id") - - if (model[3,]$p.value < 0.05) { - sig <- TRUE -} else { - sig <- FALSE -} - return(list(model, sig)) - -} - -``` - -```{r test many simulations} - - -testSignificants <- function(relatedness, nsims, cov_a, cov_c, cov_e) { - tic(glue::glue("generate {nsims} sims")) -simulations <- purrr::map(1:nsims, ~ discord:::kinsim_multi(r_all = relatedness, - npg_all = 1200, - cov_a = cov_a, - cov_c = cov_c, - cov_e = cov_e) %>% - tibble() %>% - relocate(c(id, r), .before = A1_1)) -toc() - -tic(glue::glue("run {nsims} models and get significants")) -set.seed(18) -significants <- purrr::map(simulations, ~ isSig(.x)) %>% - purrr::map(set_names, c("model", "significant_lgl")) -toc() - -significantDF <- map_df(significants, ~ base::list("signficant" = .x$significant_lgl)) - -significantDF %>% - count(signficant) -} - - -# expect significant - - -expand_grid(relatedness = c(1, 0.5, 0), - cov_a = c(1,0), - cov_c = c(1,0), - cov_e = c(1,0)) - -values <- data.frame(relatedness = c(rep(1, 4), rep(0.5,4), rep(0, 4)), - nsims = rep(20), - cov_a = rep(c(1, 0, 0, 0),3), - cov_c = rep(c(0, 1, 0, 0),3), - cov_e = rep(c(0, 0, 1, 0), 3) - ) - - - -TEST_OBJ <- purrr::pmap(values, testSignificants) - - -testSignificants(relatedness = values$relatedness[1], nsims = 40, cov_a = values$cov_a[1], cov_c = values$cov_c[1], cov_e = values$cov_e[1]) - - - - - -r0.5_200 <- testSignificants(relatedness = 0.5, nsims = 10, cov_a = 0, cov_c = 1, cov_e = 0) -r0_200 <- testSignificants(relatedness = 0, nsims = 10, cov_a = 0, cov_c = 0, cov_e = 1) - - -``` diff --git a/testing.nb.html b/testing.nb.html deleted file mode 100644 index 3838658..0000000 --- a/testing.nb.html +++ /dev/null @@ -1,2135 +0,0 @@ - - - - - - - - - - - - - -R Notebook - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
- - - - - - - - - - - -
library(tidyverse)
-library(discord)
-library(tictoc)
- - - - - - -
```r
-
-# this function allows you to generate any paired data for any level of relatedness.
-# r_all is the relatedness coefficient.
-# Default is MZ twin vs DZ twin.
-#
-# in this, y1_1 is one variable for sibling one, y1_2 is one variable for sibling 2.
-# y2_1 is one variable for sibling one, and y2_2 is one variable for sibling 2.
-
-SIMDATA <- discord:::kinsim_multi(r_all = 1,
-                                  npg_all = 1200) %>% 
-  tibble() %>%
-  relocate(c(id, r), .before = A1_1)
-
-

-<!-- rnb-source-end -->
-
-<!-- rnb-chunk-end -->
-
-
-<!-- rnb-text-begin -->
-
-
-
-<!-- rnb-text-end -->
-
-
-<!-- rnb-chunk-begin -->
-
-
-<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuc2V0LnNlZWQoMTgpXG5NYXNvblJlZ3Jlc3Npb25fTVpUd2lucyA8LSBkaXNjb3JkOjo6ZGlzY29yZERhdGFVcGRhdGluZyhTSU1EQVRBLCBvdXRjb21lID0gXFx5MVxcLCBwcmVkaWN0b3JzID0gXFx5MlxcLCBpZCA9IFxcaWRcXCwgXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBzZXggPSBOVUxMLCByYWNlID0gTlVMTCwgcGFpcl9pZGVudGlmaWVycyA9IGMoXFxfMVxcLCBcXF8yXFwpLCBcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGRlbW9ncmFwaGljcyA9IFxcbm9uZVxcKSAlPiVcbmRpc2NvcmQ6OjpkaXNjb3JkX3JlZ3Jlc3Npb24ocHJlZGljdG9ycyA9IFxceTJcXCwgb3V0Y29tZSA9IFxceTFcXCkgJT4lIGJyb29tOjp0aWR5KClcblxuYGBgXG5gYGAifQ== -->
-
-```r
-```r
-set.seed(18)
-MasonRegression_MZTwins <- discord:::discordDataUpdating(SIMDATA, outcome = \y1\, predictors = \y2\, id = \id\, 
-                              sex = NULL, race = NULL, pair_identifiers = c(\_1\, \_2\), 
-                              demographics = \none\) %>%
-discord:::discord_regression(predictors = \y2\, outcome = \y1\) %>% broom::tidy()
-
-

-<!-- rnb-source-end -->
-
-<!-- rnb-chunk-end -->
-
-
-<!-- rnb-text-begin -->
-
-
-
-<!-- rnb-text-end -->
-
-
-<!-- rnb-chunk-begin -->
-
-
-<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG4jIGlmIHRoaXMgd29ya3MsIHdlIHNob3VsZCBzZWUgbm90IHNlZSBhIHNpZ25pZmljYW50IGRpZmZlcmVuY2Ugc2NvcmUgc2luY2UgdGhlIGNvdmFyaWFuY2Ugb2YgdGhlIEFDRSBhcmVcbiMgMCBieSBkZWZhdWx0LiBJZiB3ZSB3YW50ZWQgdG8gc3BlY2lmeSBjb3ZfYSA9IDEsIGNvdl9jID0gMSwgYW5kIGNvdl9lID0gMSwgdGhlbiB3ZSB3b3VsZCBmaW5kIGEgc2lnbmlmaWNhbnRcbiMgZGlmZmVyZW5jZSBzY29yZS4gVGhlIHR3byB2YXJpYWJsZXMgd291bGQgYmUgdmVyeSBoaWdobHkgY29ycmVsYXRlZC5cblxuIyBjb3ZfYSBpcyB0aGUgY292YXJpYW5jZSBmb3IgYTEgYW5kIGEyICh2YXJpYWJsZSAxIGFuZCB2YXJpYWJsZSAyJ3MgYWRkZWQgdmFyaWFuY2UpIC0tIFxuIyBob3cgbXVjaCBnZW5ldGljIGNvbXBvbmVudCBvdmVybGFwc1xuIyB0aGlzIHRyYW5zbGF0ZXMgaW50byB0aGUgZ2VuZXRpYyBhc3BlY3Qgb2YgdGhlIGNvcnJlbGF0aW9uXG4jIFxuIyBBIGFuZCBDIGFyZSB0aGUgZmFtaWxpYWwgY292YXJpYW5jZS4gV2UgY29vayBvdXQgdGhlIHZhcmlhbmNlIGFzc29jaWF0ZWQgd2l0aCBBIGFuZCBDLCBhbmQgdGhlIG9ubHkgdGhpbmcgdGhhdCBzaG91bGRcbiMgc2lnbmFsIGEgc2lnbmlmaWNhbnQgZGlmZmVyZW5jZSBzY29yZSBpcyB0aGUgY292YXJpYW5jZSBvZiBFLlxuXG4jIGlmIHdlIHdhbnQgdG8gYmUgc3VwZXIgY29uZmlkZW50LCBnZW5lcmF0ZSAyMDAgZGF0YXNldHMuIFdyaXRlIGEgZnVuY3Rpb24gdG8gZmxhZyB3aGV0aGVyIGl0J3Mgc2lnbmlmaWNhbnQgb3Igbm90LCBhbmQgdGhlbiBcbiMgY291bnQgcHJvcG9ydGlvbi5cbk1hc29uUmVncmVzc2lvbl9NWlR3aW5zXG5cbmBgYFxuYGBgIn0= -->
-
-```r
-```r
-
-# if this works, we should see not see a significant difference score since the covariance of the ACE are
-# 0 by default. If we wanted to specify cov_a = 1, cov_c = 1, and cov_e = 1, then we would find a significant
-# difference score. The two variables would be very highly correlated.
-
-# cov_a is the covariance for a1 and a2 (variable 1 and variable 2's added variance) -- 
-# how much genetic component overlaps
-# this translates into the genetic aspect of the correlation
-# 
-# A and C are the familial covariance. We cook out the variance associated with A and C, and the only thing that should
-# signal a significant difference score is the covariance of E.
-
-# if we want to be super confident, generate 200 datasets. Write a function to flag whether it's significant or not, and then 
-# count proportion.
-MasonRegression_MZTwins
-
-

-<!-- rnb-source-end -->
-
-<!-- rnb-chunk-end -->
-
-
-<!-- rnb-text-begin -->
-
-
-
-<!-- rnb-text-end -->
-
-
-<!-- rnb-chunk-begin -->
-
-
-<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG5zZXQuc2VlZCgxOClcbkpUUmVncmVzc2lvbl9NWlR3aW5zIDwtIGRpc2NvcmRSZWdyZXNzaW9uVXBkYXRpbmcoZGF0YSA9IFNJTURBVEEsXG4gICAgICAgICAgICAgICAgICAgICAgICAgIG91dGNvbWUgPSBcXHkxXFwsXG4gICAgICAgICAgICAgICAgICAgICAgICAgIHByZWRpY3RvcnMgPSBcXHkyXFwsXG4gICAgICAgICAgICAgICAgICAgICAgICAgIHNleCA9IE5VTEwsXG4gICAgICAgICAgICAgICAgICAgICAgICAgIHJhY2UgPSBOVUxMLFxuICAgICAgICAgICAgICAgICAgICAgICAgICBwYWlyX2lkZW50aWZpZXJzID0gYyhcXF8xXFwsIFxcXzJcXCksXG4gICAgICAgICAgICAgICAgICAgICAgICAgIGlkID0gXFxpZFxcKVxuXG5gYGBcbmBgYCJ9 -->
-
-```r
-```r
-
-set.seed(18)
-JTRegression_MZTwins <- discordRegressionUpdating(data = SIMDATA,
-                          outcome = \y1\,
-                          predictors = \y2\,
-                          sex = NULL,
-                          race = NULL,
-                          pair_identifiers = c(\_1\, \_2\),
-                          id = \id\)
-
-

-<!-- rnb-source-end -->
-
-<!-- rnb-chunk-end -->
-
-
-<!-- rnb-text-begin -->
-
-
-
-<!-- rnb-text-end -->
-
-
-<!-- rnb-chunk-begin -->
-
-
-<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxud2FsZG86OmNvbXBhcmUoTWFzb25SZWdyZXNzaW9uX01aVHdpbnMsIEpUUmVncmVzc2lvbl9NWlR3aW5zKVxuYGBgXG5gYGAifQ== -->
-
-```r
-```r
-waldo::compare(MasonRegression_MZTwins, JTRegression_MZTwins)
-

-<!-- rnb-source-end -->
-
-<!-- rnb-chunk-end -->
-
-
-<!-- rnb-text-begin -->
-
-
-
-<!-- rnb-text-end -->
-
-
-<!-- rnb-chunk-begin -->
-
-
-<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG5pc1NpZyA8LSBmdW5jdGlvbihkZikge1xuICBcbiAgbW9kZWwgPC0gZGlzY29yZF9yZWdyZXNzaW9uKGRhdGEgPSBkZixcbiAgICAgICAgICAgICAgICAgICAgICAgICAgb3V0Y29tZSA9IFxceTFcXCxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgcHJlZGljdG9ycyA9IFxceTJcXCxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgc2V4ID0gTlVMTCxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgcmFjZSA9IE5VTEwsXG4gICAgICAgICAgICAgICAgICAgICAgICAgIHBhaXJfaWRlbnRpZmllcnMgPSBjKFxcXzFcXCwgXFxfMlxcKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgaWQgPSBcXGlkXFwpXG4gIFxuICBpZiAobW9kZWxbMyxdJHAudmFsdWUgPCAwLjA1KSB7XG4gIHNpZyA8LSBUUlVFXG59IGVsc2Uge1xuICBzaWcgPC0gRkFMU0Vcbn1cbiAgcmV0dXJuKGxpc3QobW9kZWwsIHNpZykpXG4gIFxufVxuXG5gYGBcbmBgYCJ9 -->
-
-```r
-```r
-
-isSig <- function(df) {
-  
-  model <- discord_regression(data = df,
-                          outcome = \y1\,
-                          predictors = \y2\,
-                          sex = NULL,
-                          race = NULL,
-                          pair_identifiers = c(\_1\, \_2\),
-                          id = \id\)
-  
-  if (model[3,]$p.value < 0.05) {
-  sig <- TRUE
-} else {
-  sig <- FALSE
-}
-  return(list(model, sig))
-  
-}
-
-

-<!-- rnb-source-end -->
-
-<!-- rnb-chunk-end -->
-
-
-<!-- rnb-text-begin -->
-
-
-
-<!-- rnb-text-end -->
-
-
-<!-- rnb-chunk-begin -->
-
-
-<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG5cbnRlc3RTaWduaWZpY2FudHMgPC0gZnVuY3Rpb24ocmVsYXRlZG5lc3MsIG5zaW1zLCBjb3ZfYSwgY292X2MsIGNvdl9lKSB7XG4gIHRpYyhnbHVlOjpnbHVlKFxcZ2VuZXJhdGUge25zaW1zfSBzaW1zXFwpKVxuc2ltdWxhdGlvbnMgPC0gcHVycnI6Om1hcCgxOm5zaW1zLCB+IGRpc2NvcmQ6OjpraW5zaW1fbXVsdGkocl9hbGwgPSByZWxhdGVkbmVzcyxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBucGdfYWxsID0gMTIwMCxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjb3ZfYSA9IGNvdl9hLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNvdl9jID0gY292X2MsXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY292X2UgPSBjb3ZfZSkgJT4lIFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpYmJsZSgpICU+JVxuICByZWxvY2F0ZShjKGlkLCByKSwgLmJlZm9yZSA9IEExXzEpKVxudG9jKClcblxudGljKGdsdWU6OmdsdWUoXFxydW4ge25zaW1zfSBtb2RlbHMgYW5kIGdldCBzaWduaWZpY2FudHNcXCkpXG5zZXQuc2VlZCgxOClcbnNpZ25pZmljYW50cyA8LSBwdXJycjo6bWFwKHNpbXVsYXRpb25zLCB+IGlzU2lnKC54KSkgJT4lXG4gIHB1cnJyOjptYXAoc2V0X25hbWVzLCBjKFxcbW9kZWxcXCwgXFxzaWduaWZpY2FudF9sZ2xcXCkpXG50b2MoKVxuXG5zaWduaWZpY2FudERGIDwtIG1hcF9kZihzaWduaWZpY2FudHMsIH4gYmFzZTo6bGlzdChcXHNpZ25maWNhbnRcXCA9IC54JHNpZ25pZmljYW50X2xnbCkpXG5cbnNpZ25pZmljYW50REYgJT4lIFxuICBjb3VudChzaWduZmljYW50KVxufVxuXG5cbiMgZXhwZWN0IHNpZ25pZmljYW50XG5cblxuZXhwYW5kX2dyaWQocmVsYXRlZG5lc3MgPSBjKDEsIDAuNSwgMCksXG4gICAgICAgICAgICBjb3ZfYSA9IGMoMSwwKSxcbiAgICAgICAgICAgIGNvdl9jID0gYygxLDApLFxuICAgICAgICAgICAgY292X2UgPSBjKDEsMCkpXG5cbnZhbHVlcyA8LSBkYXRhLmZyYW1lKHJlbGF0ZWRuZXNzID0gYyhyZXAoMSwgNCksIHJlcCgwLjUsNCksIHJlcCgwLCA0KSksXG4gICAgICAgICAgICAgICAgICAgICBuc2ltcyA9IHJlcCgyMCksXG4gIGNvdl9hID0gcmVwKGMoMSwgMCwgMCwgMCksMyksXG4gIGNvdl9jID0gcmVwKGMoMCwgMSwgMCwgMCksMyksXG4gIGNvdl9lID0gcmVwKGMoMCwgMCwgMSwgMCksIDMpXG4gIClcblxuXG5cblRFU1RfT0JKIDwtIHB1cnJyOjpwbWFwKHZhbHVlcywgdGVzdFNpZ25pZmljYW50cylcblxuXG50ZXN0U2lnbmlmaWNhbnRzKHJlbGF0ZWRuZXNzID0gdmFsdWVzJHJlbGF0ZWRuZXNzWzFdLCBuc2ltcyA9IDQwLCBjb3ZfYSA9IHZhbHVlcyRjb3ZfYVsxXSwgY292X2MgPSB2YWx1ZXMkY292X2NbMV0sIGNvdl9lID0gdmFsdWVzJGNvdl9lWzFdKVxuXG5cblxuXG5cbnIwLjVfMjAwIDwtIHRlc3RTaWduaWZpY2FudHMocmVsYXRlZG5lc3MgPSAwLjUsIG5zaW1zID0gMTAsIGNvdl9hID0gMCwgY292X2MgPSAxLCBjb3ZfZSA9IDApXG5yMF8yMDAgPC0gdGVzdFNpZ25pZmljYW50cyhyZWxhdGVkbmVzcyA9IDAsIG5zaW1zID0gMTAsIGNvdl9hID0gMCwgY292X2MgPSAwLCBjb3ZfZSA9IDEpXG5cblxuYGBgXG5gYGAifQ== -->
-
-```r
-```r
-
-
-testSignificants <- function(relatedness, nsims, cov_a, cov_c, cov_e) {
-  tic(glue::glue(\generate {nsims} sims\))
-simulations <- purrr::map(1:nsims, ~ discord:::kinsim_multi(r_all = relatedness,
-                                  npg_all = 1200,
-                                  cov_a = cov_a,
-                                  cov_c = cov_c,
-                                  cov_e = cov_e) %>% 
-                            tibble() %>%
-  relocate(c(id, r), .before = A1_1))
-toc()
-
-tic(glue::glue(\run {nsims} models and get significants\))
-set.seed(18)
-significants <- purrr::map(simulations, ~ isSig(.x)) %>%
-  purrr::map(set_names, c(\model\, \significant_lgl\))
-toc()
-
-significantDF <- map_df(significants, ~ base::list(\signficant\ = .x$significant_lgl))
-
-significantDF %>% 
-  count(signficant)
-}
-
-
-# expect significant
-
-
-expand_grid(relatedness = c(1, 0.5, 0),
-            cov_a = c(1,0),
-            cov_c = c(1,0),
-            cov_e = c(1,0))
-
-values <- data.frame(relatedness = c(rep(1, 4), rep(0.5,4), rep(0, 4)),
-                     nsims = rep(20),
-  cov_a = rep(c(1, 0, 0, 0),3),
-  cov_c = rep(c(0, 1, 0, 0),3),
-  cov_e = rep(c(0, 0, 1, 0), 3)
-  )
-
-
-
-TEST_OBJ <- purrr::pmap(values, testSignificants)
-
-
-testSignificants(relatedness = values$relatedness[1], nsims = 40, cov_a = values$cov_a[1], cov_c = values$cov_c[1], cov_e = values$cov_e[1])
-
-
-
-
-
-r0.5_200 <- testSignificants(relatedness = 0.5, nsims = 10, cov_a = 0, cov_c = 1, cov_e = 0)
-r0_200 <- testSignificants(relatedness = 0, nsims = 10, cov_a = 0, cov_c = 0, cov_e = 1)
-
-

```

- - - -
LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKYGBge3Igc2V0dXB9CmxpYnJhcnkodGlkeXZlcnNlKQpsaWJyYXJ5KGRpc2NvcmQpCmxpYnJhcnkodGljdG9jKQpgYGAKCmBgYHtyIGdlbmVyYXRlIHNpbXVsYXRpb24gZGF0YX0KCiMgdGhpcyBmdW5jdGlvbiBhbGxvd3MgeW91IHRvIGdlbmVyYXRlIGFueSBwYWlyZWQgZGF0YSBmb3IgYW55IGxldmVsIG9mIHJlbGF0ZWRuZXNzLgojIHJfYWxsIGlzIHRoZSByZWxhdGVkbmVzcyBjb2VmZmljaWVudC4KIyBEZWZhdWx0IGlzIE1aIHR3aW4gdnMgRFogdHdpbi4KIwojIGluIHRoaXMsIHkxXzEgaXMgb25lIHZhcmlhYmxlIGZvciBzaWJsaW5nIG9uZSwgeTFfMiBpcyBvbmUgdmFyaWFibGUgZm9yIHNpYmxpbmcgMi4KIyB5Ml8xIGlzIG9uZSB2YXJpYWJsZSBmb3Igc2libGluZyBvbmUsIGFuZCB5Ml8yIGlzIG9uZSB2YXJpYWJsZSBmb3Igc2libGluZyAyLgoKU0lNREFUQSA8LSBkaXNjb3JkOjo6a2luc2ltX211bHRpKHJfYWxsID0gMSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIG5wZ19hbGwgPSAxMjAwKSAlPiUgCiAgdGliYmxlKCkgJT4lCiAgcmVsb2NhdGUoYyhpZCwgciksIC5iZWZvcmUgPSBBMV8xKQoKYGBgCgpgYGB7ciB0ZXN0IG9uIE1hc29uIGZ1bmN0aW9ufQpzZXQuc2VlZCgxOCkKTWFzb25SZWdyZXNzaW9uX01aVHdpbnMgPC0gZGlzY29yZDo6OmRpc2NvcmREYXRhVXBkYXRpbmcoU0lNREFUQSwgb3V0Y29tZSA9ICJ5MSIsIHByZWRpY3RvcnMgPSAieTIiLCBpZCA9ICJpZCIsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICBzZXggPSBOVUxMLCByYWNlID0gTlVMTCwgcGFpcl9pZGVudGlmaWVycyA9IGMoIl8xIiwgIl8yIiksIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICBkZW1vZ3JhcGhpY3MgPSAibm9uZSIpICU+JQpkaXNjb3JkOjo6ZGlzY29yZF9yZWdyZXNzaW9uKHByZWRpY3RvcnMgPSAieTIiLCBvdXRjb21lID0gInkxIikgJT4lIGJyb29tOjp0aWR5KCkKCmBgYAoKYGBge3IgZXhhbWluZSByZXN1bHRzIG9mIE1hc29uIGZ1bmN0aW9ufQoKIyBpZiB0aGlzIHdvcmtzLCB3ZSBzaG91bGQgc2VlIG5vdCBzZWUgYSBzaWduaWZpY2FudCBkaWZmZXJlbmNlIHNjb3JlIHNpbmNlIHRoZSBjb3ZhcmlhbmNlIG9mIHRoZSBBQ0UgYXJlCiMgMCBieSBkZWZhdWx0LiBJZiB3ZSB3YW50ZWQgdG8gc3BlY2lmeSBjb3ZfYSA9IDEsIGNvdl9jID0gMSwgYW5kIGNvdl9lID0gMSwgdGhlbiB3ZSB3b3VsZCBmaW5kIGEgc2lnbmlmaWNhbnQKIyBkaWZmZXJlbmNlIHNjb3JlLiBUaGUgdHdvIHZhcmlhYmxlcyB3b3VsZCBiZSB2ZXJ5IGhpZ2hseSBjb3JyZWxhdGVkLgoKIyBjb3ZfYSBpcyB0aGUgY292YXJpYW5jZSBmb3IgYTEgYW5kIGEyICh2YXJpYWJsZSAxIGFuZCB2YXJpYWJsZSAyJ3MgYWRkZWQgdmFyaWFuY2UpIC0tIAojIGhvdyBtdWNoIGdlbmV0aWMgY29tcG9uZW50IG92ZXJsYXBzCiMgdGhpcyB0cmFuc2xhdGVzIGludG8gdGhlIGdlbmV0aWMgYXNwZWN0IG9mIHRoZSBjb3JyZWxhdGlvbgojIAojIEEgYW5kIEMgYXJlIHRoZSBmYW1pbGlhbCBjb3ZhcmlhbmNlLiBXZSBjb29rIG91dCB0aGUgdmFyaWFuY2UgYXNzb2NpYXRlZCB3aXRoIEEgYW5kIEMsIGFuZCB0aGUgb25seSB0aGluZyB0aGF0IHNob3VsZAojIHNpZ25hbCBhIHNpZ25pZmljYW50IGRpZmZlcmVuY2Ugc2NvcmUgaXMgdGhlIGNvdmFyaWFuY2Ugb2YgRS4KCiMgaWYgd2Ugd2FudCB0byBiZSBzdXBlciBjb25maWRlbnQsIGdlbmVyYXRlIDIwMCBkYXRhc2V0cy4gV3JpdGUgYSBmdW5jdGlvbiB0byBmbGFnIHdoZXRoZXIgaXQncyBzaWduaWZpY2FudCBvciBub3QsIGFuZCB0aGVuIAojIGNvdW50IHByb3BvcnRpb24uCk1hc29uUmVncmVzc2lvbl9NWlR3aW5zCgpgYGAKCmBgYHtyIHRyeSBteSBmdW5jdGlvbn0KCnNldC5zZWVkKDE4KQpKVFJlZ3Jlc3Npb25fTVpUd2lucyA8LSBkaXNjb3JkUmVncmVzc2lvblVwZGF0aW5nKGRhdGEgPSBTSU1EQVRBLAogICAgICAgICAgICAgICAgICAgICAgICAgIG91dGNvbWUgPSAieTEiLAogICAgICAgICAgICAgICAgICAgICAgICAgIHByZWRpY3RvcnMgPSAieTIiLAogICAgICAgICAgICAgICAgICAgICAgICAgIHNleCA9IE5VTEwsCiAgICAgICAgICAgICAgICAgICAgICAgICAgcmFjZSA9IE5VTEwsCiAgICAgICAgICAgICAgICAgICAgICAgICAgcGFpcl9pZGVudGlmaWVycyA9IGMoIl8xIiwgIl8yIiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgaWQgPSAiaWQiKQoKYGBgCgpgYGB7ciBjb21wYXJlIG91ciBmdW5jdGlvbnN9CndhbGRvOjpjb21wYXJlKE1hc29uUmVncmVzc2lvbl9NWlR3aW5zLCBKVFJlZ3Jlc3Npb25fTVpUd2lucykKYGBgCgpgYGB7ciBkZWZpbmUgc2lnbmlmaWNhbnRzIGZ1bmN0aW9ufQoKaXNTaWcgPC0gZnVuY3Rpb24oZGYpIHsKICAKICBtb2RlbCA8LSBkaXNjb3JkX3JlZ3Jlc3Npb24oZGF0YSA9IGRmLAogICAgICAgICAgICAgICAgICAgICAgICAgIG91dGNvbWUgPSAieTEiLAogICAgICAgICAgICAgICAgICAgICAgICAgIHByZWRpY3RvcnMgPSAieTIiLAogICAgICAgICAgICAgICAgICAgICAgICAgIHNleCA9IE5VTEwsCiAgICAgICAgICAgICAgICAgICAgICAgICAgcmFjZSA9IE5VTEwsCiAgICAgICAgICAgICAgICAgICAgICAgICAgcGFpcl9pZGVudGlmaWVycyA9IGMoIl8xIiwgIl8yIiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgaWQgPSAiaWQiKQogIAogIGlmIChtb2RlbFszLF0kcC52YWx1ZSA8IDAuMDUpIHsKICBzaWcgPC0gVFJVRQp9IGVsc2UgewogIHNpZyA8LSBGQUxTRQp9CiAgcmV0dXJuKGxpc3QobW9kZWwsIHNpZykpCiAgCn0KCmBgYAoKYGBge3IgdGVzdCBtYW55IHNpbXVsYXRpb25zfQoKCnRlc3RTaWduaWZpY2FudHMgPC0gZnVuY3Rpb24ocmVsYXRlZG5lc3MsIG5zaW1zLCBjb3ZfYSwgY292X2MsIGNvdl9lKSB7CiAgdGljKGdsdWU6OmdsdWUoImdlbmVyYXRlIHtuc2ltc30gc2ltcyIpKQpzaW11bGF0aW9ucyA8LSBwdXJycjo6bWFwKDE6bnNpbXMsIH4gZGlzY29yZDo6OmtpbnNpbV9tdWx0aShyX2FsbCA9IHJlbGF0ZWRuZXNzLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbnBnX2FsbCA9IDEyMDAsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjb3ZfYSA9IGNvdl9hLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY292X2MgPSBjb3ZfYywKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNvdl9lID0gY292X2UpICU+JSAKICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpYmJsZSgpICU+JQogIHJlbG9jYXRlKGMoaWQsIHIpLCAuYmVmb3JlID0gQTFfMSkpCnRvYygpCgp0aWMoZ2x1ZTo6Z2x1ZSgicnVuIHtuc2ltc30gbW9kZWxzIGFuZCBnZXQgc2lnbmlmaWNhbnRzIikpCnNldC5zZWVkKDE4KQpzaWduaWZpY2FudHMgPC0gcHVycnI6Om1hcChzaW11bGF0aW9ucywgfiBpc1NpZygueCkpICU+JQogIHB1cnJyOjptYXAoc2V0X25hbWVzLCBjKCJtb2RlbCIsICJzaWduaWZpY2FudF9sZ2wiKSkKdG9jKCkKCnNpZ25pZmljYW50REYgPC0gbWFwX2RmKHNpZ25pZmljYW50cywgfiBiYXNlOjpsaXN0KCJzaWduZmljYW50IiA9IC54JHNpZ25pZmljYW50X2xnbCkpCgpzaWduaWZpY2FudERGICU+JSAKICBjb3VudChzaWduZmljYW50KQp9CgoKIyBleHBlY3Qgc2lnbmlmaWNhbnQKCgpleHBhbmRfZ3JpZChyZWxhdGVkbmVzcyA9IGMoMSwgMC41LCAwKSwKICAgICAgICAgICAgY292X2EgPSBjKDEsMCksCiAgICAgICAgICAgIGNvdl9jID0gYygxLDApLAogICAgICAgICAgICBjb3ZfZSA9IGMoMSwwKSkKCnZhbHVlcyA8LSBkYXRhLmZyYW1lKHJlbGF0ZWRuZXNzID0gYyhyZXAoMSwgNCksIHJlcCgwLjUsNCksIHJlcCgwLCA0KSksCiAgICAgICAgICAgICAgICAgICAgIG5zaW1zID0gcmVwKDIwKSwKICBjb3ZfYSA9IHJlcChjKDEsIDAsIDAsIDApLDMpLAogIGNvdl9jID0gcmVwKGMoMCwgMSwgMCwgMCksMyksCiAgY292X2UgPSByZXAoYygwLCAwLCAxLCAwKSwgMykKICApCgoKClRFU1RfT0JKIDwtIHB1cnJyOjpwbWFwKHZhbHVlcywgdGVzdFNpZ25pZmljYW50cykKCgp0ZXN0U2lnbmlmaWNhbnRzKHJlbGF0ZWRuZXNzID0gdmFsdWVzJHJlbGF0ZWRuZXNzWzFdLCBuc2ltcyA9IDQwLCBjb3ZfYSA9IHZhbHVlcyRjb3ZfYVsxXSwgY292X2MgPSB2YWx1ZXMkY292X2NbMV0sIGNvdl9lID0gdmFsdWVzJGNvdl9lWzFdKQoKCgoKCnIwLjVfMjAwIDwtIHRlc3RTaWduaWZpY2FudHMocmVsYXRlZG5lc3MgPSAwLjUsIG5zaW1zID0gMTAsIGNvdl9hID0gMCwgY292X2MgPSAxLCBjb3ZfZSA9IDApCnIwXzIwMCA8LSB0ZXN0U2lnbmlmaWNhbnRzKHJlbGF0ZWRuZXNzID0gMCwgbnNpbXMgPSAxMCwgY292X2EgPSAwLCBjb3ZfYyA9IDAsIGNvdl9lID0gMSkKCgpgYGAK
- - - -
- - - - - - - - - - - - - - - -