11 | | Example code |
| 11 | Example code importing only counts pseudotime and cells weights on each trajectory |
| 12 | |
| 13 | {{{ |
| 14 | #! /usr/bin/Rscript |
| 15 | |
| 16 | library(tradeSeq) |
| 17 | library(monocle) |
| 18 | library(RColorBrewer) |
| 19 | library(SingleCellExperiment) |
| 20 | library(slingshot) |
| 21 | RNGversion("3.5.0") |
| 22 | palette(brewer.pal(8, "Dark2")) |
| 23 | |
| 24 | packageDescription("tradeSeq")$Version |
| 25 | #load data |
| 26 | |
| 27 | cellPSMatrix <- as.matrix(read.table("Subset.PsTime.txt", sep="\t")) |
| 28 | cellWeights <- as.matrix(read.table("Subset.CellsWeight.txt", sep="\t")) |
| 29 | counts <- as.matrix(read.table("RawCounts.txt", sep="\t")) |
| 30 | |
| 31 | #this step is optional |
| 32 | #select the subset of cells that are present on the pseudotime and cells weight file |
| 33 | countsFiltered <- subset.matrix(counts, select=rownames(cellWeights)) |
| 34 | |
| 35 | #Evaluate the knots that we'll use when fitting the General Additive Models |
| 36 | icMat <- evaluateK(counts = countsFiltered, k = 3:10, |
| 37 | nGenes = 200, verbose = T, pseudotime = cellPSMatrix, |
| 38 | cellWeights = cellWeights) |
| 39 | |
| 40 | library(BiocParallel) |
| 41 | library(magrittr) |
| 42 | |
| 43 | # Register BiocParallel Serial Execution (no parallelization in that case) |
| 44 | BiocParallel::register(BiocParallel::SerialParam()) |
| 45 | FitModels <- fitGAM(counts = countsFiltered,nknots = 8, pseudotime = cellPSMatrix, cellWeights = cellWeights ) |
| 46 | |
| 47 | patternRes <- patternTest(FitModels) |
| 48 | oPat <- order(patternRes$waldStat, decreasing = TRUE) |
| 49 | head(rownames(patternRes)[oPat]) |
| 50 | save.image("Data_afterPatternTest.RData") |
| 51 | |
| 52 | pdf( "plotSmoothers.pdf", w=11, h=8.5) |
| 53 | plotSmoothers(FitModels, countsFiltered , gene = rownames(patternRes)[oPat][1]) |
| 54 | plotSmoothers(FitModels, countsFiltered , gene = rownames(patternRes)[oPat][2]) |
| 55 | plotSmoothers(FitModels, countsFiltered , gene = rownames(patternRes)[oPat][3]) |
| 56 | plotSmoothers(FitModels, countsFiltered , gene = rownames(patternRes)[oPat][4]) |
| 57 | plotSmoothers(FitModels, countsFiltered , gene = rownames(patternRes)[oPat][5]) |
| 58 | dev.off() |
| 59 | |
| 60 | write.table(patternRes, file="PatternTest8knots.txt", sep="\t", quote=F, row.names=T) |
| 61 | |
| 62 | |
| 63 | }}} |
| 64 | |