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

Factor in amplification efficiency #21

Open richardstoeckl opened 4 years ago

richardstoeckl commented 4 years ago

Hi, thank you for your very useful package!

In your publication you mention the ability to factor in less-than-perfect amplification efficiencies in a future version of this package. Sadly, I was not able to find this option. Is this feature still in development?

Greetings, Richard

richardstoeckl commented 4 years ago

So I tried a workaround myself, maybe you might have an Idea how to implement it better.

I modified the "pcr_ddct" function to also take a dataframe of amplification efficiencies as Input (which I reformatted from the output of pcr_assess to only include the genes as columns with one row where the amplification efficiency is stored).

I then added a variable "y" to the function that sets "res", which uses the variable "x" to get the column name of the column with the values in "x":


res <- apply(goi,
               MARGIN = 2,
               FUN = function(x) {
                 # get the name of the column that x points to
                 y = names(goi)[which(goi == x, arr.ind=T)[, "col"]]
                 # calculate the averages
                 x_ave <- .pcr_average(x, group_var)
                 ref_ave <- .pcr_average(ref, group_var)

and called a modified version of the ".pcr_relative" function with the addition of the amplification_efficiency dataframe and the column name:


                 # calculate the relative expression
                 rel_expr <- .pcr_relative_mod(ddct, amplification_efficiency, y)

                 # calculate the error bars
                 upper <- .pcr_relative_mod(ddct - error, amplification_efficiency, y)
                 lower <- .pcr_relative_mod(ddct + error, amplification_efficiency, y)

The modified function looks like this:


.pcr_relative_mod <- function(vec, amp_eff, y) {
  a_e <- amp_eff %>% select(contains(y))
  res <- a_e[1,1] ^ (-vec)
  return(res)
}

It uses the column name (gene name) stored in "y" to get the associated amplification efficiency of the gene from the dataframe and uses that as the base instead of the "2" in the original function.

What are your thoughts on this approach?

Best regards, Richard

MahShaaban commented 3 years ago

Thanks, @richardstoeckl for reaching out, and sorry for the late reply. I think this approach would work. Would you like to submit a pull request? I will review the code and make sure it works with existing functions.