Loading packages

Let’s start by leading the packages needed for the analysis. We will need to tidyverse, that get loaded in the following code chunk.

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.2.1     ✔ readr     2.2.0
## ✔ forcats   1.0.1     ✔ stringr   1.6.0
## ✔ ggplot2   4.0.3     ✔ tibble    3.3.1
## ✔ lubridate 1.9.5     ✔ tidyr     1.3.2
## ✔ purrr     1.2.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

Loading data

Let’s load the mouse RNA sequencing data from Blackmore et al. (2017), The effect of upper-respiratory infection on transcriptomic changes in the CNS.

rna <- read.csv("data/rnaseq.csv")
rna
##      gene     sample expression     organism age    sex  infection  strain time
## 1     Asl GSM2545336       1170 Mus musculus   8 Female InfluenzaA C57BL/6    8
## 2    Apod GSM2545336      36194 Mus musculus   8 Female InfluenzaA C57BL/6    8
## 3 Cyp2d22 GSM2545336       4060 Mus musculus   8 Female InfluenzaA C57BL/6    8
## 4    Klk6 GSM2545336        287 Mus musculus   8 Female InfluenzaA C57BL/6    8
## 5   Fcrls GSM2545336         85 Mus musculus   8 Female InfluenzaA C57BL/6    8
##       tissue mouse ENTREZID
## 1 Cerebellum    14   109900
## 2 Cerebellum    14    11815
## 3 Cerebellum    14    56448
## 4 Cerebellum    14    19144
## 5 Cerebellum    14    80891
##                                                                        product
## 1                               argininosuccinate lyase, transcript variant X1
## 2                                       apolipoprotein D, transcript variant 3
## 3 cytochrome P450, family 2, subfamily d, polypeptide 22, transcript variant 2
## 4                         kallikrein related-peptidase 6, transcript variant 2
## 5                Fc receptor-like S, scavenger receptor, transcript variant X1
##      ensembl_gene_id external_synonym chromosome_name   gene_biotype
## 1 ENSMUSG00000025533    2510006M18Rik               5 protein_coding
## 2 ENSMUSG00000022548             <NA>              16 protein_coding
## 3 ENSMUSG00000061740             2D22              15 protein_coding
## 4 ENSMUSG00000050063             Bssp               7 protein_coding
## 5 ENSMUSG00000015852    2810439C17Rik               3 protein_coding
##                             phenotype_description
## 1           abnormal circulating amino acid level
## 2                      abnormal lipid homeostasis
## 3                        abnormal skin morphology
## 4                         abnormal cytokine level
## 5 decreased CD8-positive alpha-beta T cell number
##   hsapiens_homolog_associated_gene_name
## 1                                   ASL
## 2                                  APOD
## 3                                CYP2D6
## 4                                  KLK6
## 5                                 FCRL2
##  [ reached 'max' / getOption("max.print") -- omitted 32423 rows ]

Visualisation

We love gglot2, don’t we?

rna |> 
  filter(gene %in% c("Dok3", "Asl", "Pank4", 
                      "Trim72", "Gabrr2")) |> 
  group_by(gene, time, sex) |> 
  summarise(meanexp = mean(expression)) |> 
  ggplot(aes(x = time, y = meanexp, colour = sex)) +
  geom_line() +
  facet_wrap(~ gene)
## `summarise()` has regrouped the output.
## ℹ Summaries were computed grouped by gene, time, and sex.
## ℹ Output is grouped by gene and time.
## ℹ Use `summarise(.groups = "drop_last")` to silence this message.
## ℹ Use `summarise(.by = c(gene, time, sex))` for per-operation grouping
##   (`?dplyr::dplyr_by`) instead.

Session information

sessionInfo()
## R version 4.6.1 (2026-06-24)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.4 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so;  LAPACK version 3.12.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: Etc/UTC
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] lubridate_1.9.5 forcats_1.0.1   stringr_1.6.0   dplyr_1.2.1    
##  [5] purrr_1.2.2     readr_2.2.0     tidyr_1.3.2     tibble_3.3.1   
##  [9] ggplot2_4.0.3   tidyverse_2.0.0
## 
## loaded via a namespace (and not attached):
##  [1] gtable_0.3.6       jsonlite_2.0.0     compiler_4.6.1     tidyselect_1.2.1  
##  [5] jquerylib_0.1.4    scales_1.4.0       yaml_2.3.12        fastmap_1.2.0     
##  [9] R6_2.6.1           labeling_0.4.3     generics_0.1.4     knitr_1.51        
## [13] bslib_0.11.0       pillar_1.11.1      RColorBrewer_1.1-3 tzdb_0.5.0        
## [17] rlang_1.3.0        stringi_1.8.7      cachem_1.1.0       xfun_0.60         
## [21] sass_0.4.10        S7_0.2.2           otel_0.2.0         timechange_0.4.0  
## [25] cli_3.6.6          withr_3.0.3        magrittr_2.0.5     digest_0.6.39     
## [29] grid_4.6.1         rstudioapi_0.19.0  hms_1.1.4          lifecycle_1.0.5   
## [33] vctrs_0.7.3        evaluate_1.0.5     glue_1.8.1         farver_2.1.2      
## [37] rmarkdown_2.31     tools_4.6.1        pkgconfig_2.0.3    htmltools_0.5.9