Construct heatmap for non-steady state (NSS) analysis with improved mechanism score

This uses the output of bakR and a differential expression analysis software to construct a dataframe that can be passed to pheatmap::pheatmap(). This heatmap will display the result of a steady-state quasi-independent analysis of NR-seq data.

Usage

NSSHeat2(
  bakRFit,
  DE_df,
  bakRModel = c("MLE", "Hybrid", "MCMC"),
  DE_cutoff = 0.05,
  bakR_cutoff = 0.05,
  Exp_ID = 2,
  a = 0.01
)

Arguments

bakRFit: bakRFit object
DE_df: dataframe of required format with differential expression analysis results. See Further-Analyses vignette for details on what this dataframe should look like
bakRModel: Model fit from which bakR implementation should be used? Options are MLE, Hybrid, or MCMC
DE_cutoff: padj cutoff for calling a gene differentially expressed
bakR_cutoff: padj cutoff for calling a fraction new significantly changed
Exp_ID: Exp_ID of experimental sample whose comparison to the reference sample you want to use. Only one reference vs. experimental sample comparison can be used at a time
a: Scale factor for mechanism scoring function. Larger numbers will make it so increase in mechanism score as a function of differential expression or differential fraction new z-score is more dramatic.

Value

returns data frame that can be passed to pheatmap::pheatmap()

Details

Unlike NSSHeat, NSSHeat2 uses a mechanism scoring function that is not discontinuous as the "degradation driven" vs. "synthesis driven" boundary. Instead, the score approaches 0 as the function approaches the boundary from either side.

Examples

# \donttest{
# Simulate small dataset
sim <- Simulate_bakRData(100, nreps = 2)

# Analyze data with bakRFit
Fit <- bakRFit(sim$bakRData)
#> Finding reliable Features
#> Filtering out unwanted or unreliable features
#> Processing data...
#> Estimating labeled mutation rate
#> Estimated pnews for each sample are:
#> # A tibble: 4 × 3
#> # Groups:   mut [2]
#>     mut  reps   pnew
#>   <int> <dbl>  <dbl>
#> 1     1     1 0.0498
#> 2     1     2 0.0497
#> 3     2     1 0.0498
#> 4     2     2 0.0497
#> Estimating unlabeled mutation rate
#> Estimated pold is:  0.001
#> Estimating fraction labeled
#> Estimating per replicate uncertainties
#> Estimating read count-variance relationship
#> Averaging replicate data and regularizing estimates
#> Assessing statistical significance
#> All done! Run QC_checks() on your bakRFit object to assess the quality of your data and get recommendations for next steps.

# Number of features that made it past filtering
NF <- nrow(Fit$Fast_Fit$Effects_df)

# Simulate mock differential expression data frame
DE_df <- data.frame(XF = as.character(1:NF),
                       log2FoldChange = stats::rnorm(NF, 0, 2))

DE_df$stat <- DE_df$log2FoldChange/0.5

DE_df$padj <- 2*stats::dnorm(-abs(DE_df$stat))

# make heatmap matrix
Heatmap <- NSSHeat2(bakRFit = Fit,
               DE_df = DE_df,
               bakRModel = "MLE")

# }