Browse code

added option to shrink F stat denominator using limma::squeezeVar

Nils Kurzawa authored on 09/04/2020 08:24:10
Showing5 changed files

... ...
@@ -1,7 +1,7 @@
1 1
 Package: TPP2D
2 2
 Title: Detection of ligand-protein interactions from 2D thermal 
3 3
     profiles (DLPTP)
4
-Version: 1.3.13
4
+Version: 1.3.14
5 5
 Authors@R: c(
6 6
     person("Nils", "Kurzawa", email = "nils.kurzawa@embl.de", 
7 7
         role = c("aut", "cre")),
... ...
@@ -23,5 +23,5 @@ URL: http://bioconductor.org/packages/TPP2D
23 23
 RoxygenNote: 6.1.1
24 24
 Depends: R (>= 3.6.0), stats, utils, dplyr, methods
25 25
 Imports: ggplot2, tidyr, foreach, doParallel, openxlsx, stringr, 
26
-  RCurl, parallel, MASS, BiocParallel
26
+  RCurl, parallel, MASS, BiocParallel, limma
27 27
 Suggests: knitr, testthat
... ...
@@ -37,6 +37,7 @@ importFrom(RCurl,url.exists)
37 37
 importFrom(doParallel,registerDoParallel)
38 38
 importFrom(foreach,"%dopar%")
39 39
 importFrom(foreach,foreach)
40
+importFrom(limma,squeezeVar)
40 41
 importFrom(methods,is)
41 42
 importFrom(methods,new)
42 43
 importFrom(openxlsx,read.xlsx)
... ...
@@ -5,6 +5,8 @@
5 5
 #' fitAndEvalDataset
6 6
 #' @param df_null data frame containing results from analysis by
7 7
 #' bootstrapNull
8
+#' @param squeezeDenominator logical indicating whether F statistic
9
+#' denominator should be shrinked using limma::squeezeVar
8 10
 #' 
9 11
 #' @return data frame annotating each protein with a FDR based on 
10 12
 #' it's F statistic and number of observations
... ...
@@ -23,12 +25,17 @@
23 25
 #' @export
24 26
 #'
25 27
 #' @import dplyr
26
-getFDR <- function(df_out, df_null){
28
+getFDR <- function(df_out, df_null, squeezeDenominator = TRUE){
27 29
     dataset <- nObs <- nObsRound <- F_statistic <- 
28 30
         is_decoy <- max_rank <- true_cumsum <- 
29 31
         null_cumsum <- representative <- clustername <- 
30 32
         dataset <- FDR <- all_true <- all_null <- NULL
31 33
     
34
+    if(squeezeDenominator){
35
+        df_out <- .shrinkFstat(df_out, trueOrNull = "true")
36
+        df_null <- .shrinkFstat(df_null, trueOrNull = "null")
37
+    }
38
+    
32 39
     B <- max(as.numeric(
33 40
         gsub("bootstrap_", "", unique(df_null$dataset))))
34 41
     
... ...
@@ -52,6 +59,31 @@ getFDR <- function(df_out, df_null){
52 59
     return(out_df)
53 60
 }
54 61
 
62
+#' @importFrom limma squeezeVar 
63
+.shrinkFstat <- function(inDf, trueOrNull = "true"){
64
+    rssH1 <- df2 <- rssH0 <- rssH1Squeezed <- df1 <- 
65
+        nObs <- nObsRound <- dataset <- NULL
66
+    if(trueOrNull == "true"){
67
+        outDf <- inDf %>% 
68
+            mutate(nObsRound = round(nObs/10)*10) %>% 
69
+            group_by(nObsRound) %>% 
70
+            mutate(rssH1Squeezed = limma::squeezeVar(
71
+                rssH1, df = df2)$var.post) %>% 
72
+            ungroup %>% 
73
+            mutate(F_statistic = (rssH0 - rssH1)/
74
+                       (rssH1Squeezed) * df2/df1)
75
+    }else if(trueOrNull == "null"){
76
+        outDf <- inDf %>% 
77
+            mutate(nObsRound = round(nObs/10)*10) %>% 
78
+            group_by(dataset, nObsRound) %>% 
79
+            mutate(rssH1Squeezed = limma::squeezeVar(
80
+                rssH1, df = df2)$var.post) %>% 
81
+            ungroup %>% 
82
+            mutate(F_statistic = (rssH0 - rssH1)/
83
+                       (rssH1Squeezed) * df2/df1)
84
+    }
85
+    return(outDf)
86
+}
55 87
 #' Compute FDR for given F statistics based on true and
56 88
 #' null dataset (old function)
57 89
 #' 
... ...
@@ -5,7 +5,7 @@
5 5
 \title{Get FDR for given F statistics based on true and
6 6
 null dataset}
7 7
 \usage{
8
-getFDR(df_out, df_null)
8
+getFDR(df_out, df_null, squeezeDenominator = TRUE)
9 9
 }
10 10
 \arguments{
11 11
 \item{df_out}{data frame containing results from analysis by
... ...
@@ -13,6 +13,9 @@ fitAndEvalDataset}
13 13
 
14 14
 \item{df_null}{data frame containing results from analysis by
15 15
 bootstrapNull}
16
+
17
+\item{squeezeDenominator}{logical indicating whether F statistic
18
+denominator should be shrinked using limma::squeezeVar}
16 19
 }
17 20
 \value{
18 21
 data frame annotating each protein with a FDR based on 
... ...
@@ -7,4 +7,4 @@ test_that("computeFStatFromParams works as expected", {
7 7
                                       filter(representative %in% 1:3))
8 8
     fstat_df <- computeFStatFromParams(params_df)
9 9
     expect_equal(round(fstat_df$F_statistic, 2), c(1.23, 0.36, 0.24), tolerance = 0.015)
10
-})
11 10
\ No newline at end of file
11
+})