1 . Overall study design

We collected kidneys from 188 Diversity Outbred mice (93 females and 95 males) at 6, 12, and 18 months of age (63, 62, and 63 animals respectively). The right kidney was homogenized and an aliquot was used for RNA-sequencing and another for shotgun proteomics. A week prior to sacrifice, we also collected urine samples for renal physiological phenotyping (albumin, phosphate, and creatinine levels).

2 . Data access

We have provided a copy of our data and code on our GitHub repository.

3 . Changes in mRNA and protein with age

We fitted a linear model for all 6,667 genes with both mRNA and protein expression information.

\[ y_{i} \sim Age_{i} + Sex_{i} + Generation_{i} \] From linear model, we recorded slope and p-value generated from each mRNA and protein in the downstream processes. Additionally for a subset of genes with both mRNA and protein information, we also extracted the slope coefficient from.

# Make sure correct paths and used and paste below into cluster to run:

qsub -v I="./Rdata/DO188b_kidney_noprobs.RData ./Anova_output/kidney_anova_slope_output.csv",script=anova_tests_slope_pairs Rsubmit_args.sh

4 . Direction of Change with Age

We can illustrate the direction of change in mRNA and protein expression with age using the slope calculated from the linear model in the previous section.

In our manuscript, we used clusterProfiler to identify enriched biological functions and pathways of each of the 4 significant cluster of genes. For more information please refer to our manuscript HERE. (needs link)

5 . The Ecological Fallacy

Previously we showed how mRNA and protein expression changed with age. We want to highlight a few genes to explain what is known as “The ecological fallacy”, whereby the correlation of the total population is not the same as individual age group correlation.

We picked Slc5a12 and Flot1 as the best examples to demonstrate this ecological fallacy.


At the whole population level, the mRNA expression of Slc5a12 decreased with age and the protein expression increased with age. However, when broken down by age groups, the mRNA and protein expressions are positively correlated. When were compare between the age groups we do not see a shift between 6 month and 12 months, but at 18 months of age the average of the group shifts to the top right. The upward shift indicates that the protein expression increased but the leftward shift indicates the decrease of mRNA expression.

In contrast to Slc5a12, Flot1 is an example where by the mRNA expression increases with age and the protein expression decreases with age at the whole population level. At the level of the individual age groups, we observe a bottom-right shift with an increase in age. The downward shift indicates a decrease in protein expression, and a rightward shift indicates a increase in mRNA expression with age.

6 . Age-related changes in protein expression are not mediated by mRNA expression

In order to further clarify the role of transcriptional regulation in determining age-related changes in protein, we carried a mediation analysis between protein and its corresponding mRNA for 6,667 genes.

We found that age-related changes in protein levels are not mediated by mRNA expression, however counter to this, many sex-related changes in protein levels are regulated by mRNA expression.

7 . Relating to renal physiology

In addition to gene expression and protein expression quantification, we also measured metabolites in the urine samples (albumin, phosphate, and creatinine) collected from the same 188 animals 1 week prior to their sacrifice. From the analysis conducted in previous sections, we identified Slc34a1, Slc34a3, and Lrp2 in group D to be of particular interest. Genes in group D have decreasing mRNA expression and increasing protein expression increasing with age at the whole population level.

8 .Drivers of transcriptomic and proteomic changes.

Using mRNA expression (eQTL) and protein expression (pQTL) as phenotypes in a quantitative trait loci (QTL) mapping approach, we were able to identify regions that regulate age- and sex-associated changes in gene expression.
For more information on QTL analyses, Karl Broman’s R/qtl2 user guide is great resource for understanding and getting started with QTL scans.

# Load QTL2 library
library(qtl2)
library(qtl2convert)

# Load probs data
load("./RNAseq_data/DO188b_kidney.RData")

# Reformat from DOQTL format to meet qtl2 format specification
probs <- probs_doqtl_to_qtl2(genoprobs, snps, pos_column = "bp")
snps$chr <- as.character(snps$chr)
snps$chr[snps$chr=="X"] <- "20" #change to numeric for creating map later
map <- map_df_to_list(map = snps, pos_column = "bp")

# Create matrix for each model
## Additive
addcovar <- model.matrix(~ Sex + Age + Generation, data=annot.samples)
## Sex-interactive
sex_intcovar <- model.matrix(~ Sex, data = annot.samples)
## Age-interactive
age_intcovar <- model.matrix(~ Age, data = annot.samples)

# scan additive model:
  LODadd <- scan1(genoprobs=probs,
               kinship=Glist,
               pheno=expr.protein[,p],
               addcovar=addcovar[,-1],
               cores=3,
               reml=TRUE)

# scan full model for age-interaction:
  LODint <- scan1(genoprobs=probs,
               kinship=Glist,
               pheno=expr.protein[,p],
               addcovar=addcovar[,-1],
               intcovar=age_intcovar[,-1], # replace with sex_intcovar for sex-interaction
               cores=3,
               reml=TRUE)

# get marker with maximum LOD score difference between full and additive scans:
  maxMarker <- function(full, add){
    diff <- as.data.frame(full)
    colnames(diff) <- "FullLOD"
    diff$AddLOD <- add[,1]
    diff$IntAgeLODDiff <- diff$FullLOD - diff$AddLOD
    max <- diff[which(diff$IntAgeLODDiff == max(diff$IntAgeLODDiff, na.rm = TRUE)[1])[1],]
    max$IntAgeChr <- str_split_fixed(rownames(max),"_",2)[,1]
    max$IntAgePos <- snps[snps$marker == rownames(max),]$bp
    return(max)
  }
  Max <- maxMarker(full = LODint, add = LODadd)
}

# output
output <- output %>% mutate(
  FullLOD = Max$FullLOD,
  AddLOD = Max$AddLOD,
  IntAgeChr = Max$IntAgeChr,
  IntAgePos = Max$IntAgePos,
  IntAgeLODDiff = Max$IntAgeLODDiff
)

9 . Session Information

## R version 3.5.1 (2018-07-02)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: macOS  10.14
## 
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] grid      stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] scales_1.0.0    bindrcpp_0.2.2  ggsci_2.9       gridExtra_2.3  
##  [5] forcats_0.3.0   stringr_1.3.1   dplyr_0.7.6     purrr_0.2.5    
##  [9] readr_1.1.1     tidyr_0.8.1     tibble_1.4.2    ggplot2_3.0.0  
## [13] tidyverse_1.2.1
## 
## loaded via a namespace (and not attached):
##  [1] Rcpp_0.12.19     cellranger_1.1.0 pillar_1.3.0     compiler_3.5.1  
##  [5] plyr_1.8.4       bindr_0.1.1      tools_3.5.1      digest_0.6.17   
##  [9] lubridate_1.7.4  jsonlite_1.5     evaluate_0.11    nlme_3.1-137    
## [13] gtable_0.2.0     lattice_0.20-35  pkgconfig_2.0.2  rlang_0.2.2     
## [17] cli_1.0.1        rstudioapi_0.8   yaml_2.2.0       haven_1.1.2     
## [21] withr_2.1.2      xml2_1.2.0       httr_1.3.1       knitr_1.20      
## [25] hms_0.4.2        rprojroot_1.3-2  tidyselect_0.2.4 glue_1.3.0      
## [29] R6_2.3.0         readxl_1.1.0     rmarkdown_1.10   modelr_0.1.2    
## [33] magrittr_1.5     backports_1.1.2  htmltools_0.3.6  rvest_0.3.2     
## [37] assertthat_0.2.0 colorspace_1.3-2 labeling_0.3     stringi_1.2.4   
## [41] lazyeval_0.2.1   munsell_0.5.0    broom_0.5.0      crayon_1.3.4

Last updated: 25 October, 2018