MahShaaban / pcr

Quality assessing, analyzing and testing the statistical significance of real-time quantitative PCR data
https://CRAN.R-project.org/package=pcr
GNU General Public License v3.0
27 stars 8 forks source link

pcr_test function no longer provides estimate using "wilcox.test" #3

Closed whatever347 closed 6 years ago

whatever347 commented 6 years ago

Since a "recent" (a few months I guess) update, the pcr_test function from the package no longer provides an estimate in the resulting data.table using the "wilcox.test" argument. The "estimate" column now only contains NA.

Reproducible example using the vignette:

fl <- system.file('extdata', 'ct4.csv', package = 'pcr')
ct4 <- readr::read_csv(fl)

group <- rep(c('control', 'treatment'), each = 12)

Results <- pcr_test(ct4,
         group_var = group,
         reference_gene = 'ref',
         reference_group = 'control',
         test = 'wilcox.test')

Results

I made a workaround by manually extracting a estimate from the result of the call to the base R wilcox.test(..., conf.int = T) function, but this quickly gets tedious. I quickly parsed the code of the pcr_test and the pcr_wilcox function i guess that the problems arises from this part of the pcr_wilcox function which doesn't properly extract the estimate value from the htest object produced by pcr_wilcox :

wilcox_test <- data_frame(estimate = unname(wilcox_test$estimate[1] - 
            wilcox_test$estimate[2]), p_value = wilcox_test$p.value, 
            lower = wilcox_test$conf.int[1], upper = wilcox_test$conf.int[2])

Thanks for this otherwise great package.

MahShaaban commented 6 years ago

Thanks for @whatever347 for catching this. You are right. wilcox_test should be as follows:

wilcox_test <- data_frame(
      estimate = unname(wilcox_test$estimate), # the base wilcox.test retruns a single estimate
      p_value = wilcox_test$p.value,
      lower = wilcox_test$conf.int[1],
      upper = wilcox_test$conf.int[2]
                        )

Using the reproducible example you provided, I get the expected results.

> fl <- system.file('extdata', 'ct4.csv', package = 'pcr')
> ct4 <- readr::read_csv(fl)
Parsed with column specification:
cols(
  ref = col_double(),
  target = col_double()
)
> 
> group <- rep(c('control', 'treatment'), each = 12)
> 
> Results <- pcr_test(ct4,
+                     group_var = group,
+                     reference_gene = 'ref',
+                     reference_group = 'control',
+                     test = 'wilcox.test')
> 
> Results
# A tibble: 1 x 5
  gene   estimate    p_value  lower  upper
  <chr>     <dbl>      <dbl>  <dbl>  <dbl>
1 target   -0.635 0.00000296 -0.881 -0.423

I made the necessary change. You can use the GitHub version of pcr for now, and the next CRAN release will include the fix. To install the package from GitHub, use

devtools::install_github('MahShaaban/pcr')

Please let me know if you have any further comments.