Description of your vignette
Data comprises 10 .RCC files of 5 subjects, 2 Timepoint per subjects
Note that the values of the QC metrics vary only slightly from sample to sample. This is one indication of good data.
Get the normalised counts
# install.packages("pheatmap", "RColorBrewer", "viridis")
library(pheatmap)
library(RColorBrewer)
library(viridis)
# Data frame with column annotations.
mat_col <- data.frame(group = col_groups)
rownames(mat_col) <- colnames(mat)
# List with colors for each annotation.
mat_colors <- list(group = brewer.pal(3, "Set1"))
names(mat_colors$group) <- unique(col_groups)
pheatmap(
mat = mat,
color = inferno(10),
border_color = NA,
show_colnames = FALSE,
show_rownames = FALSE,
annotation_col = mat_col,
annotation_colors = mat_colors,
drop_levels = TRUE,
fontsize = 14,
main = "Default Heatmap"
)
FOV registration QC Each individual lane scanned on an nCounter system is divided into a few hundred imaging sections, called Fields of View (FOV). The system images these FOVs separately, and sums the barcode counts (representing the GeoMx DSP probe counts) of all FOVs. Finally, the system reports the number of FOVs successfully imaged as FOV Counted. Significant discrepancy between the number of FOV for which imaging was attempted (FOV Count) and for which imaging was successful (FOV Counted) may indicate an issue with imaging performance. Recommended percentage of registered FOVs (i.e., FOV Counted over FOV Count) is 75%. Lanes will be flagged if this percentage is lower.
p<- ggsummarystats(
df, x = "gene", y = "Count_Norm",
ggfunc = ggboxplot, add = c("jitter"),
color = "Timepoint", palette = "npg"
)
p %>% ggplotly()
Make sure count matrix and phenotypes have the same samples
samples_kept <- intersect(selected_pheno[["sampleid"]], colnames(expr_counts))
expr_counts <- expr_counts[, samples_kept]
selected_pheno <- filter(selected_pheno, sampleid %in% !!samples_kept)
selected_pheno$Timepoint <- factor(selected_pheno$Timepoint, levels = unique(selected_pheno$Timepoint))
selected_pheno$Subject <- factor(selected_pheno$Subject, levels = unique(selected_pheno$Subject))
design <- model.matrix(~ Subject + Timepoint , selected_pheno)
eBayes(lmFit(expr_counts, design))
## Warning: Zero sample variances detected, have been offset away from zero
## An object of class "MArrayLM"
## $coefficients
## (Intercept) Subject2 Subject3 Subject4 Subject5
## TNFSF13B 302.5 -89.5 186.5 -46.5 -57.0
## CCL5 77.0 -49.5 -43.0 -47.5 -37.0
## MCL1 3727.5 171.5 1801.0 719.5 -243.5
## PPBP 145.3 -104.5 -102.0 -89.5 -106.5
## CCL20 24.9 -12.5 -1.0 0.5 -18.5
## TimepointPost-Treatment
## TNFSF13B -48.0
## CCL5 -24.0
## MCL1 -552.0
## PPBP -46.6
## CCL20 -12.8
## 125 more rows ...
##
## $rank
## [1] 6
##
## $assign
## [1] 0 1 1 1 1 2
##
## $qr
## $qr
## (Intercept) Subject2 Subject3 Subject4 Subject5
## 1 -3.1622777 -0.6324555 -0.63245553 -0.6324555 -0.63245553
## 2 0.3162278 1.2649111 -0.31622777 -0.3162278 -0.31622777
## 3 0.3162278 -0.6704429 1.22474487 -0.4082483 -0.40824829
## 4 0.3162278 -0.6704429 0.04378207 1.1547005 -0.57735027
## 5 0.3162278 0.1201265 -0.67804553 -0.0526541 -1.00000000
## 6 0.3162278 0.1201265 -0.67804553 -0.0526541 0.08609256
## 7 0.3162278 0.1201265 0.13845105 -0.6784390 0.10928778
## 8 0.3162278 0.1201265 0.13845105 -0.6784390 0.10928778
## 9 0.3162278 0.1201265 0.13845105 0.1875864 0.69328512
## 10 0.3162278 0.1201265 0.13845105 0.1875864 0.69328512
## TimepointPost-Treatment
## 1 -1.581139e+00
## 2 1.110223e-16
## 3 8.326673e-17
## 4 3.330669e-16
## 5 -2.775558e-16
## 6 -1.581139e+00
## 7 1.457241e-01
## 8 7.781797e-01
## 9 -2.300550e-01
## 10 4.024005e-01
## attr(,"assign")
## [1] 0 1 1 1 1 2
## attr(,"contrasts")
## attr(,"contrasts")$Subject
## [1] "contr.treatment"
##
## attr(,"contrasts")$Timepoint
## [1] "contr.treatment"
##
##
## $qraux
## [1] 1.316228 1.120127 1.043782 1.059320 1.086093 1.397931
##
## $pivot
## [1] 1 2 3 4 5 6
##
## $tol
## [1] 1e-07
##
## $rank
## [1] 6
##
##
## $df.residual
## [1] 4 4 4 4 4
## 125 more elements ...
##
## $sigma
## TNFSF13B CCL5 MCL1 PPBP CCL20
## 153.866013 14.361407 1704.431650 36.798777 8.879752
## 125 more elements ...
##
## $cov.coefficients
## (Intercept) Subject2 Subject3 Subject4
## (Intercept) 0.6 -5.000000e-01 -5.000000e-01 -5.000000e-01
## Subject2 -0.5 1.000000e+00 5.000000e-01 5.000000e-01
## Subject3 -0.5 5.000000e-01 1.000000e+00 5.000000e-01
## Subject4 -0.5 5.000000e-01 5.000000e-01 1.000000e+00
## Subject5 -0.5 5.000000e-01 5.000000e-01 5.000000e-01
## TimepointPost-Treatment -0.2 -1.358774e-16 -1.211652e-16 -1.708889e-16
## Subject5 TimepointPost-Treatment
## (Intercept) -5.000000e-01 -2.000000e-01
## Subject2 5.000000e-01 -1.358774e-16
## Subject3 5.000000e-01 -1.211652e-16
## Subject4 5.000000e-01 -1.708889e-16
## Subject5 1.000000e+00 -1.110223e-16
## TimepointPost-Treatment -1.110223e-16 4.000000e-01
##
## $stdev.unscaled
## (Intercept) Subject2 Subject3 Subject4 Subject5
## TNFSF13B 0.7745967 1 1 1 1
## CCL5 0.7745967 1 1 1 1
## MCL1 0.7745967 1 1 1 1
## PPBP 0.7745967 1 1 1 1
## CCL20 0.7745967 1 1 1 1
## TimepointPost-Treatment
## TNFSF13B 0.6324555
## CCL5 0.6324555
## MCL1 0.6324555
## PPBP 0.6324555
## CCL20 0.6324555
## 125 more rows ...
##
## $pivot
## [1] 1 2 3 4 5 6
##
## $Amean
## TNFSF13B CCL5 MCL1 PPBP CCL20
## 277.2 29.6 3941.2 41.5 12.2
## 125 more elements ...
##
## $method
## [1] "ls"
##
## $design
## (Intercept) Subject2 Subject3 Subject4 Subject5 TimepointPost-Treatment
## 1 1 0 0 0 0 0
## 2 1 0 0 0 0 1
## 3 1 1 0 0 0 0
## 4 1 1 0 0 0 1
## 5 1 0 1 0 0 0
## 6 1 0 1 0 0 1
## 7 1 0 0 1 0 0
## 8 1 0 0 1 0 1
## 9 1 0 0 0 1 0
## 10 1 0 0 0 1 1
## attr(,"assign")
## [1] 0 1 1 1 1 2
## attr(,"contrasts")
## attr(,"contrasts")$Subject
## [1] "contr.treatment"
##
## attr(,"contrasts")$Timepoint
## [1] "contr.treatment"
##
##
## $df.prior
## [1] 0.5057102
##
## $s2.prior
## [1] 305.6438
##
## $var.prior
## [1] 5.234852e-02 3.271783e-05 3.271783e-05 3.271783e-05 3.271783e-05
## [6] 5.234852e-02
##
## $proportion
## [1] 0.01
##
## $s2.post
## TNFSF13B CCL5 MCL1 PPBP CCL20
## 21051.8572 217.4057 2579061.4834 1236.4682 104.3048
## 125 more elements ...
##
## $t
## (Intercept) Subject2 Subject3 Subject4 Subject5
## TNFSF13B 2.691563 -0.6168475 1.28538619 -0.32048503 -0.3928526
## CCL5 6.741858 -3.3571420 -2.91630519 -3.22149992 -2.5093789
## MCL1 2.996479 0.1067907 1.12145762 0.44802263 -0.1516241
## PPBP 5.334559 -2.9718359 -2.90073933 -2.54525657 -3.0287131
## CCL20 3.147542 -1.2239337 -0.09791469 0.04895735 -1.8114218
## TimepointPost-Treatment
## TNFSF13B -0.5230775
## CCL5 -2.5736279
## MCL1 -0.5434733
## PPBP -2.0953880
## CCL20 -1.9816540
## 125 more rows ...
##
## $df.total
## [1] 4.50571 4.50571 4.50571 4.50571 4.50571
## 125 more elements ...
##
## $p.value
## (Intercept) Subject2 Subject3 Subject4 Subject5
## TNFSF13B 0.048136660 0.56715434 0.26075275 0.76292284 0.71231726
## CCL5 0.001626842 0.02365487 0.03756737 0.02718176 0.05919489
## MCL1 0.034455412 0.91954427 0.31827659 0.67484756 0.88603478
## PPBP 0.004225330 0.03537952 0.03820800 0.05680989 0.03328790
## CCL20 0.029358250 0.28109480 0.92620417 0.96304811 0.13622056
## TimepointPost-Treatment
## TNFSF13B 0.62560274
## CCL5 0.05499968
## MCL1 0.61258922
## PPBP 0.09639318
## CCL20 0.11061633
## 125 more rows ...
##
## $lods
## (Intercept) Subject2 Subject3 Subject4 Subject5
## TNFSF13B -4.497262 -4.595129 -4.595112 -4.595134 -4.595133
## CCL5 -4.428247 -4.595072 -4.595077 -4.595073 -4.595084
## MCL1 -4.485776 -4.595136 -4.595117 -4.595132 -4.595136
## PPBP -4.439306 -4.595077 -4.595078 -4.595083 -4.595076
## CCL20 -4.480749 -4.595114 -4.595136 -4.595136 -4.595098
## TimepointPost-Treatment
## TNFSF13B -4.638315
## CCL5 -4.460168
## MCL1 -4.636945
## PPBP -4.494717
## CCL20 -4.504111
## 125 more rows ...
##
## $F
## [1] 6.892100 10.363879 10.429779 5.284191 3.995662
## 125 more elements ...
##
## $F.p.value
## [1] 0.03227335 0.01451789 0.01433408 0.05284917 0.08625346
## 125 more elements ...
library(purrr)
res <- MYDF[["nacho"]] %>%
# filter( QC PARAMETER) %>% # possible additional QC filters
filter(grepl("Endogenous", CodeClass)) %>%
dplyr::select(sampleid, Name, Accession, Count, Count_Norm, Timepoint) %>%
mutate_all(~ na_if(.x, "NA")) %>%
drop_na() %>%
group_by(Name, Accession) %>%
nest() %>%
ungroup() %>%
slice(1:10) %>% # the ten first genes
mutate(
lm = map(.x = data, .f = function(idata) {# lm or whatever model you want
as.data.frame(summary(lm(formula = Count_Norm ~ Timepoint, idata))$coef)
})
)%>%
dplyr::select(-data)
unnest(res, lm)
mutate(
lm = map(.x = data, .f = function(idata) {# lm or whatever model you want
as.data.frame(summary(lm(formula = Count_Norm ~ `disease.event:ch1`, idata))$coef)
})
) %>%
select(-data)
## R version 4.1.2 (2021-11-01)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 22000)
##
## Matrix products: default
##
## locale:
## [1] LC_COLLATE=English_United States.1252
## [2] LC_CTYPE=English_United States.1252
## [3] LC_MONETARY=English_United States.1252
## [4] LC_NUMERIC=C
## [5] LC_TIME=English_United States.1252
##
## attached base packages:
## [1] grid stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] rstatix_0.7.0 ggpubr_0.4.0
## [3] heatmaply_1.3.0 viridis_0.6.2
## [5] viridisLite_0.4.0 plotly_4.10.0
## [7] nanostringr_0.2.0 tweeDEseq_1.40.0
## [9] tweeDEseqCountData_1.32.0 devtools_2.4.3
## [11] usethis_2.1.5 VennDiagram_1.7.1
## [13] futile.logger_1.4.3 mixOmics_6.18.1
## [15] lattice_0.20-45 MASS_7.3-54
## [17] NACHO_1.1.0 GGally_2.1.2
## [19] reshape2_1.4.4 gridExtra_2.3
## [21] RColorBrewer_1.1-2 org.Mm.eg.db_3.14.0
## [23] AnnotationDbi_1.56.2 gplots_3.1.1
## [25] Glimma_2.4.0 DESeq2_1.34.0
## [27] SummarizedExperiment_1.24.0 MatrixGenerics_1.6.0
## [29] matrixStats_0.61.0 GenomicRanges_1.46.1
## [31] GenomeInfoDb_1.30.1 IRanges_2.28.0
## [33] edgeR_3.36.0 limma_3.50.1
## [35] ggthemes_4.2.4 forcats_0.5.1
## [37] stringr_1.4.0 dplyr_1.0.8
## [39] purrr_0.3.4 readr_2.1.2
## [41] tidyr_1.1.4 tibble_3.1.6
## [43] tidyverse_1.3.1 openxlsx_4.2.5
## [45] tictoc_1.0.1 readxl_1.3.1
## [47] NanoStringNCTools_1.2.0 ggplot2_3.3.5
## [49] S4Vectors_0.32.3 NanoStringQCPro_1.26.0
## [51] Biobase_2.54.0 BiocGenerics_0.40.0
## [53] BiocStyle_2.22.0
##
## loaded via a namespace (and not attached):
## [1] pkgmaker_0.32.2 bit64_4.0.5 knitr_1.37
## [4] DelayedArray_0.20.0 data.table_1.14.2 KEGGREST_1.34.0
## [7] RCurl_1.98-1.6 doParallel_1.0.17 generics_0.1.2
## [10] cowplot_1.1.1 callr_3.7.0 lambda.r_1.2.4
## [13] RSQLite_2.2.10 bit_4.0.4 tzdb_0.2.0
## [16] webshot_0.5.2 xml2_1.3.3 lubridate_1.8.0
## [19] httpuv_1.6.5 ggsci_2.9 assertthat_0.2.1
## [22] xfun_0.29 hms_1.1.1 jquerylib_0.1.4
## [25] evaluate_0.15 promises_1.2.0.1 TSP_1.2-0
## [28] fansi_0.5.0 caTools_1.18.2 dendextend_1.15.2
## [31] dbplyr_2.1.1 igraph_1.2.11 DBI_1.1.2
## [34] geneplotter_1.72.0 htmlwidgets_1.5.4 reshape_0.8.8
## [37] rARPACK_0.11-0 ellipsis_0.3.2 RSpectra_0.16-0
## [40] crosstalk_1.2.0 backports_1.4.1 bookdown_0.24
## [43] ggiraph_0.8.2 annotate_1.72.0 gridBase_0.4-7
## [46] vctrs_0.3.8 remotes_2.4.2 abind_1.4-5
## [49] cachem_1.0.6 withr_2.4.3 prettyunits_1.1.1
## [52] cluster_2.1.2 pacman_0.5.1 lazyeval_0.2.2
## [55] crayon_1.5.0 ellipse_0.4.2 genefilter_1.76.0
## [58] pkgconfig_2.0.3 labeling_0.4.2 vipor_0.4.5
## [61] pkgload_1.2.4 seriation_1.3.2 rlang_1.0.1
## [64] lifecycle_1.0.1 registry_0.5-1 modelr_0.1.8
## [67] cellranger_1.1.0 rprojroot_2.0.2 rngtools_1.5.2
## [70] Matrix_1.3-4 carData_3.0-5 reprex_2.0.1
## [73] beeswarm_0.4.0 processx_3.5.2 pheatmap_1.0.12
## [76] png_0.1-7 bitops_1.0-7 KernSmooth_2.23-20
## [79] Biostrings_2.62.0 blob_1.2.2 nor1mix_1.3-0
## [82] ggsignif_0.6.3 scales_1.1.1 memoise_2.0.1
## [85] magrittr_2.0.1 plyr_1.8.6 zlibbioc_1.40.0
## [88] compiler_4.1.2 cli_3.1.0 XVector_0.34.0
## [91] ps_1.6.0 formatR_1.11 tidyselect_1.1.2
## [94] stringi_1.7.6 highr_0.9 yaml_2.2.2
## [97] locfit_1.5-9.4 ggrepel_0.9.1 sass_0.4.0
## [100] tools_4.1.2 parallel_4.1.2 rstudioapi_0.13
## [103] uuid_1.0-3 foreach_1.5.2 farver_2.1.0
## [106] digest_0.6.29 BiocManager_1.30.16 shiny_1.7.1
## [109] Rcpp_1.0.7 car_3.0-12 broom_0.7.12
## [112] later_1.3.0 shinyWidgets_0.6.4 org.Hs.eg.db_3.14.0
## [115] httr_1.4.2 colorspace_2.0-2 rvest_1.0.2
## [118] brio_1.1.3 XML_3.99-0.9 fs_1.5.2
## [121] splines_4.1.2 sessioninfo_1.2.2 systemfonts_1.0.3
## [124] xtable_1.8-4 jsonlite_1.7.3 futile.options_1.0.1
## [127] corpcor_1.6.10 testthat_3.1.2 R6_2.5.1
## [130] pillar_1.7.0 htmltools_0.5.2 mime_0.12
## [133] NMF_0.23.0 glue_1.5.1 fastmap_1.1.0
## [136] DT_0.21 BiocParallel_1.28.3 codetools_0.2-18
## [139] pkgbuild_1.3.1 utf8_1.2.2 bslib_0.3.1
## [142] ggbeeswarm_0.6.0 gtools_3.9.2 magick_2.7.3
## [145] zip_2.2.0 survival_3.2-13 rmarkdown_2.12
## [148] desc_1.4.0 munsell_0.5.0 GenomeInfoDbData_1.2.7
## [151] iterators_1.0.14 haven_2.4.3 gtable_0.3.0
## [154] cqn_1.40.0