## ----eval=FALSE---------------------------------------------------------------
#  # Install from github
#  library(devtools)
#  install_github("brennanhilton/RNAseqCovarImpute")
#  
#  # Install from Bioconductor (not yet on Bioconductor)
#  
#  if (!require("BiocManager", quietly = TRUE)) {
#      install.packages("BiocManager")
#  }
#  
#  BiocManager::install("RNAseqCovarImpute")

## -----------------------------------------------------------------------------
library(RNAseqCovarImpute)
library(dplyr)
library(BiocParallel)
data(example_data)
data(example_DGE)

## -----------------------------------------------------------------------------
intervals <- get_gene_bin_intervals(example_DGE, example_data, n = 10)

## ----message=FALSE, warning=FALSE, echo=FALSE---------------------------------
intervals %>%
    head(10) %>%
    knitr::kable(digits = 3, caption = "The first 10 gene bins. Start and end columns indicate row numbers for the beginning and end of each bin. Number indicates the number of genes in each bin.")

## -----------------------------------------------------------------------------
# Randomize the order of gene identifiers
annot <- example_DGE$genes
annot <- annot[sample(seq_len(nrow(annot))), ]
# Match order of the genes in the DGE to the randomized order of genes in the annotation
example_DGE <- example_DGE[annot, ]

## -----------------------------------------------------------------------------
gene_bin_impute <- impute_by_gene_bin(example_data,
    intervals,
    example_DGE,
    m = 3
)

## -----------------------------------------------------------------------------
gene_bin_impute <- impute_by_gene_bin(example_data,
    intervals,
    example_DGE,
    m = 3,
    BPPARAM = SerialParam()
)

## -----------------------------------------------------------------------------
coef_se <- limmavoom_imputed_data_list(
    gene_intervals = intervals,
    DGE = example_DGE,
    imputed_data_list = gene_bin_impute,
    m = 3,
    voom_formula = "~x + y + z + a + b"
)

## -----------------------------------------------------------------------------
final_res <- combine_rubins(
    DGE = example_DGE,
    model_results = coef_se,
    predictor = "x"
)

## ----message=FALSE, warning=FALSE, echo=FALSE---------------------------------
final_res %>%
    dplyr::select(probe, coef_combined, combined_p_bayes, combined_p_adj_bayes) %>%
    head(10) %>%
    knitr::kable(digits = 3, caption = "The top 10 genes associated with predictor x sorted by lowest P-value")

## ----sessionInfo, echo=FALSE--------------------------------------------------
sessionInfo()