inst/examples/development.r

In mariev/chromatoplotsgui: Interactive graphics and interface for chromatoplots

library(cranvas)
library(chromatoplots)
library(faahKO)
library(xcms)
setwd('/home/marie/chromatoplotsgui/R')
#source('qplotRaw.r')
#source('qplotChrom.r')

# example from chromatoplots
cdfpath <- system.file('cdf',package='faahKO')
cdffiles <- list.files(cdfpath,recursive=TRUE,full.names=TRUE)
raw1 <- loadSample(cdffiles[1])
raw_prof <- genProfile(raw1)
cor_prof <- removeBaseline(raw_prof, "median", scanrad = 100)
peaks <- findPeaks(cor_prof, "gauss")        


## re-installation
remove.packages('chromatoplotsgui')


setwd("/Users/marie/Documents/repos")


install.packages('chromatoplotsgui', repos = NULL, type = 'source', 
                 INSTALL_opts = '--no-multiarch')

## example1: using gui
library(chromatoplotsgui)
library(chromatoplots)
# setwd('/home/marie/chromatoplotsgui/data/raw')
 setwd("/Users/marie/Documents/repos/chromatoplotsgui/data/raw")
setwd('treatment_0')
raw1 <- loadSample(dir())
raw2 <- xcmsRaw(dir())
raw3 <- xcmsRaw(dir(), profparam = list(baselevel = 150, basespace = 0))
prof <- genProfile(raw1)

cpGUI()
cpGUI(raw1)
cpGUI(raw2)
cpGUI(raw3)
# cpgtest <- chromatoplotsgui(dir())
# setwd('/home/marie/chromatoplotsgui/inst/examples')
# setwd("/Users/marie/Documents/repos/chromatoplotsgui/data")
# save(cpgtest, file = 'testdata.rda')


load('/Users/marie/Documents/repos/chromatoplotsgui/data/testdata.rda')
chromatoplotsgui(guidata = cpgtest)

#### example 2: independent plots
library(chromatoplotsgui)
library(chromatoplots)
setwd("/Users/marie/Documents/repos/chromatoplotsgui/data/raw")
setwd('treatment_0')
raw1 <- loadSample(dir())
prof <- genProfile(raw1)
cor <- removeBaseline(prof, 'median', scanrad = 100)
# cor <- removeBaseline(prof, "rbe", b = 4) # document the b value
# peaks <- findPeaks(cor, 'gauss')
# save(peaks, file = 'peaks.rda')
load("/Users/marie/Documents/repos/chromatoplotsgui/data/peaks.rda")
# setwd('../treatment_1')
# ci_exp <- loadExperiment(dir()[1:2], pipeline = pipeline(peaks))
# save(ci_exp, file = 'comps.rda')
load("/Users/marie/Documents/repos/chromatoplotsgui/data/comps.rda")
xset_comps <- findComps(ci_exp, "sigma")


cgRawPlot(raw1)
cgRawPlot(xcmsRaw(dir()))

cgProfPlot(prof)
cgProfPlot(xcmsRaw(dir()))

cgRmBasePlot(cor, 112)
raw2 <- xcmsRaw(dir(), profmethod = "binlinbase", 
                profparam = list(baselevel = 150, basespace = 0))
prof2 <- raw2@env$profile
profMedFilt(raw2, 0, 100)
raw2@env$profile <- prof2 - raw2@env$profile ## calculate residuals
cgRmBasePlot(raw2, 201)

cgfindPeaksPlot(peaks, cor, 201)
cgfindPeaksPlot(peaks, cor, 209) ## zoom on leftmost peaks
cgfindPeaksPlot(peaks, cor, 110)

prof3 <- raw2@env$profile
quant <- quantile(prof3, 0.99, na.rm = TRUE)
## fit_profile is from chromatoplots
fits <- apply(prof3, 1, fit_profile, raw2@scantime, quant, TRUE)
fits_flat <- unlist(fits, F)
mz <- rep(profMz(raw2), sapply(fits, length))
good_fits <- !sapply(fits_flat, is.character)
fits <- fits_flat[good_fits]
names(fits) <- NULL
coeff_names <- c("mu", "h", "sigma", "t")
integ <- function(xy, a, b) sum((xy[(a + 1):b, 2] + xy[a:(b - 
                1), 2]) * (xy[(a + 1):b, 1] - xy[a:(b - 1), 1]))/2
peaks1 <- t(sapply(fits, function(model) {
  e <- model$m$getEnv()
  span <- diff(range(e$x))
  coeffs <- coefficients(model)[coeff_names]
  coeffs[is.na(coeffs)] <- 0
  coeffs[["sigma"]] <- abs(coeffs[["sigma"]])
  names(coeffs) <- coeff_names
  c(rt = coeffs[["mu"]], rtmin = e$x[1], rtmax = tail(e$x, 
          1), into = integ(cbind(e$x, e$y), 1, length(e$x)), 
    intf = sqrt(2 * pi) * coeffs[["h"]] * coeffs[["sigma"]], 
    maxo = max(e$y), maxf = coeffs[["h"]], tau = coeffs[["t"]], 
    sigma = coeffs[["sigma"]], error = sqrt(deviance(model)/df.residual(model)), 
    span = span/5, delta_i = coeffs[["h"]] - max(e$y), 
    delta_t = coeffs[["mu"]] - e$x[which.max(e$y)], peak_sd = sd(e$y), 
    `sigma/h` = coeffs[["sigma"]]/coeffs[["h"]])
}))
masses <- mz[good_fits]
peaks2 <- cbind(mz = masses, mzmin = masses, mzmax = masses, 
                 peaks1)
cgfindPeaksPlot(peaks2, raw2, 201)