PROGENy pathway signatures: Application to Bulk transcriptomics

Michael Schubert

2020-04-27

PROGENy pathway signatures

This R package provides the model we inferred in the publication “Perturbation-response genes reveal signaling footprints in cancer gene expression” and a function to obtain pathway scores from a gene expression matrix. It is available on bioRxiv.

Scoring the airway package data for pathway scores

This is to outline how to prepare expression data, in this case from the airway package for pathway activity analysis using PROGENy.

Checking for differences between the groups

So now we might be interested how the treatment with dexamethasone affects signaling pathways. To do this, we check if the control is different to the perturbed condition using a linear model:

## # A tibble: 14 x 5
##    estimate std.error statistic   p.value pathway 
##       <dbl>     <dbl>     <dbl>     <dbl> <chr>   
##  1    8.77      0.834    10.5   0.0000433 Androgen
##  2    1.25      3.22      0.389 0.711     EGFR    
##  3   -7.05      0.835    -8.44  0.000151  Estrogen
##  4    0.783     1.92      0.408 0.698     Hypoxia 
##  5   -0.972     0.657    -1.48  0.190     JAK-STAT
##  6    1.04      1.35      0.773 0.469     MAPK    
##  7    0.456     0.697     0.655 0.537     NFkB    
##  8   -4.42      0.789    -5.61  0.00137   p53     
##  9   -0.618     1.94     -0.319 0.761     PI3K    
## 10    2.34      1.09      2.15  0.0756    TGFb    
## 11    0.715     0.681     1.05  0.335     TNFa    
## 12    0.279     0.785     0.356 0.734     Trail   
## 13    0.728     0.700     1.04  0.338     VEGF    
## 14   -1.76      0.564    -3.12  0.0206    WNT

What we see is that indeed the p53/DNA damage response pathway is less active after treatment than before.

Reproducing drug associations on the GDSC panel

Below is an example on how to calculate pathway scores for cell lines in the Genomics of Drug Sensitivity in Cancer (GDSC) panel, and to check for associations with drug response.

The code used for the analyses is available on Github.

Getting the data

This example shows how to use the GDSC gene expression data of multiple cell lines together with PROGENy to calculate pathway activity and then to check for associations with drug sensitivity.

First, we need the GDSC data for both gene expression and drug response. They are available on the GDSC1000 web site:

You can also download the files manually (adjust the file names when loading):

Running PROGENy to get pathway activity scores

Activity inference is done using a weighted sum of the model genes. We can run this without worrying about the order of genes in the expression matrix using:

To visualize the progeny result

We now have the pathway activity scores for the pathways defined in PROGENy:

##           Androgen       EGFR     Estrogen     Hypoxia    JAK-STAT        MAPK
## 906826  -0.4241254  0.1385179 -0.269674710 -0.06875973 -0.37598907  0.41134699
## 687983  -1.7393225 -0.9213671 -0.763490640 -1.32374232 -0.94662180 -0.90704066
## 910927  -1.5393075 -0.0695914  1.613512221 -0.76482836 -1.05428545  0.17412564
## 1240138  0.9107488 -0.2200467  0.219584368 -0.73615734 -0.07536995 -0.60238190
## 1240139 -1.0602951 -0.1846655 -0.556250681  0.13664228 -0.58644833 -0.04550095
## 906792   0.7518034  0.8045299  0.008914293  0.40118664 -0.50038060  0.77199885
##               NFkB        p53       PI3K       TGFb       TNFa      Trail
## 906826  -0.6374683 -1.0903096 -0.1953806 -0.7793323 -0.4329185 -0.8018626
## 687983  -1.3157127 -0.8692247  0.2542425 -0.8242281 -1.0195560 -0.5104074
## 910927  -0.3059911  0.8896616 -0.7487845  0.7580845 -0.1756034 -0.8092882
## 1240138  0.9755970  2.0632167  0.4839630  1.9108565  0.9865983 -0.2377919
## 1240139 -0.4751547  0.8763231  1.3814461  0.9999634 -0.3916422 -0.9224493
## 906792  -0.1878324 -0.8895494 -3.0379834 -0.3389734 -0.1436331  0.2455631
##                VEGF         WNT
## 906826  -0.13639016  0.01831199
## 687983  -0.08663407 -0.72986989
## 910927   0.13321227  2.28602257
## 1240138  0.01391272  1.51533198
## 1240139  0.90296847  1.71036886
## 906792   0.18149340 -1.22160214

Testing if MAPK activity is significantly associated with Trametinib

Trametinib is a MEK inhibitor, so we would assume that cell lines that have a higher MAPK activity are more sensitive to MEK inhibition.

We can test this the following way:

## 
## Call:
## lm(formula = trametinib ~ mapk)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -6.548 -1.255  0.338  1.338  6.819 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.94703    0.06442  -14.70   <2e-16 ***
## mapk        -1.35978    0.06449  -21.09   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.998 on 960 degrees of freedom
## Multiple R-squared:  0.3165, Adjusted R-squared:  0.3158 
## F-statistic: 444.6 on 1 and 960 DF,  p-value: < 2.2e-16

And indeed we find that MAPK activity is strongly associated with sensitivity to Trametinib: the Pr(>|t|) is much smaller than the conventional threshold of 0.05.

The intercept is significant as well, but we’re not really interested if the mean drug response is above or below 0 in this case.

Note, however, that we tested all cell lines at once and did not adjust for the effect different tissues may have.

R version information

## R version 4.0.0 (2020-04-24)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.4 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.11-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.11-bioc/R/lib/libRlapack.so
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=C              
##  [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       
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] BiocFileCache_1.12.0        dbplyr_1.4.3               
##  [3] biomaRt_2.44.0              DESeq2_1.28.0              
##  [5] airway_1.7.0                SummarizedExperiment_1.18.0
##  [7] DelayedArray_0.14.0         matrixStats_0.56.0         
##  [9] Biobase_2.48.0              GenomicRanges_1.40.0       
## [11] GenomeInfoDb_1.24.0         IRanges_2.22.0             
## [13] S4Vectors_0.26.0            BiocGenerics_0.34.0        
## [15] knitr_1.28                  tibble_3.0.1               
## [17] pheatmap_1.0.12             readr_1.3.1                
## [19] tidyr_1.0.2                 ggplot2_3.3.0              
## [21] Seurat_3.1.5                dplyr_0.8.5                
## [23] progeny_1.10.0              BiocStyle_2.16.0           
## 
## loaded via a namespace (and not attached):
##   [1] readxl_1.3.1           backports_1.1.6        plyr_1.8.6            
##   [4] igraph_1.2.5           lazyeval_0.2.2         splines_4.0.0         
##   [7] BiocParallel_1.22.0    listenv_0.8.0          digest_0.6.25         
##  [10] htmltools_0.4.0        magick_2.3             fansi_0.4.1           
##  [13] gdata_2.18.0           magrittr_1.5           memoise_1.1.0         
##  [16] cluster_2.1.0          ROCR_1.0-7             limma_3.44.0          
##  [19] globals_0.12.5         annotate_1.66.0        askpass_1.1           
##  [22] prettyunits_1.1.1      colorspace_1.4-1       blob_1.2.1            
##  [25] rappdirs_0.3.1         ggrepel_0.8.2          xfun_0.13             
##  [28] crayon_1.3.4           RCurl_1.98-1.2         jsonlite_1.6.1        
##  [31] genefilter_1.70.0      survival_3.1-12        zoo_1.8-7             
##  [34] ape_5.3                glue_1.4.0             gtable_0.3.0          
##  [37] zlibbioc_1.34.0        XVector_0.28.0         leiden_0.3.3          
##  [40] future.apply_1.5.0     scales_1.1.0           DBI_1.1.0             
##  [43] Rcpp_1.0.4.6           viridisLite_0.3.0      xtable_1.8-4          
##  [46] progress_1.2.2         reticulate_1.15        bit_1.1-15.2          
##  [49] rsvd_1.0.3             tsne_0.1-3             htmlwidgets_1.5.1     
##  [52] httr_1.4.1             gplots_3.0.3           RColorBrewer_1.1-2    
##  [55] ellipsis_0.3.0         ica_1.0-2              pkgconfig_2.0.3       
##  [58] XML_3.99-0.3           farver_2.0.3           uwot_0.1.8            
##  [61] utf8_1.1.4             locfit_1.5-9.4         tidyselect_1.0.0      
##  [64] labeling_0.3           rlang_0.4.5            reshape2_1.4.4        
##  [67] AnnotationDbi_1.50.0   cellranger_1.1.0       munsell_0.5.0         
##  [70] tools_4.0.0            cli_2.0.2              generics_0.0.2        
##  [73] RSQLite_2.2.0          broom_0.5.6            ggridges_0.5.2        
##  [76] evaluate_0.14          stringr_1.4.0          yaml_2.2.1            
##  [79] npsurv_0.4-0           bit64_0.9-7            fitdistrplus_1.0-14   
##  [82] caTools_1.18.0         purrr_0.3.4            RANN_2.6.1            
##  [85] pbapply_1.4-2          future_1.17.0          nlme_3.1-147          
##  [88] compiler_4.0.0         curl_4.3               plotly_4.9.2.1        
##  [91] png_0.1-7              lsei_1.2-0             geneplotter_1.66.0    
##  [94] stringi_1.4.6          RSpectra_0.16-0        lattice_0.20-41       
##  [97] Matrix_1.2-18          vctrs_0.2.4            pillar_1.4.3          
## [100] lifecycle_0.2.0        BiocManager_1.30.10    lmtest_0.9-37         
## [103] RcppAnnoy_0.0.16       data.table_1.12.8      cowplot_1.0.0         
## [106] bitops_1.0-6           irlba_2.3.3            patchwork_1.0.0       
## [109] R6_2.4.1               bookdown_0.18          KernSmooth_2.23-17    
## [112] gridExtra_2.3          codetools_0.2-16       MASS_7.3-51.6         
## [115] gtools_3.8.2           assertthat_0.2.1       openssl_1.4.1         
## [118] withr_2.2.0            sctransform_0.2.1      GenomeInfoDbData_1.2.3
## [121] hms_0.5.3              grid_4.0.0             rmarkdown_2.1         
## [124] Rtsne_0.15