load data and annotation
dcount = readRDS('count.rds')
ano = readRDS('ano.rds')
groups = readRDS('fraction.rds')
p21 = basicP2proc(dcount,min.cells.per.gene = 0,n.cores = 5)
84 cells, 53735 genes; normalizing ...
Using plain model
Winsorizing ...
log scale ...
done.
calculating variance fit ...
using gam
101 overdispersed genes ... 101
persisting ...
done.
running PCA using 3000 OD genes .
.
Warning in irlba(x, nv = nPcs, nu = 0, center = cm, right_only = FALSE, :
You're computing too large a percentage of total singular values, use a standard svd instead.
.
.
done
creating space of type angular done
adding data ... done
building index ... done
querying ... done
Estimating embeddings.
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
*****************************************************|
*****perplexity is too large, reducing to 27
running tSNE using 5 cores:
Warning in irlba(x, nv = nPcs, nu = 0, center = cm, right_only = FALSE, :
You're computing too large a percentage of total singular values, use a standard svd instead.
creating space of type angular done
adding data ... done
building index ... done
querying ... done
df = read.csv('CellCycleGeneList_1134.txt',sep='\t',header=T)
#gsano = readRDS('gene.ano.rds')
alls = NULL
gl = list()
for (i in c("M/G1" , "G2" , "G2/M" , "G1/S" , "S phase")) {
tmp1 = df[df$PHASE == i, ]
gs = apply(tmp1, 1, function(x)
strsplit(x['NAME'], ' ')[[1]][1])
gs = gs[gs != '']
index12 = (gs %in% gsano$Gene_Name)
gs = gs[index12]
index12 = match(gs, gsano$Gene_Name)
gs2 = Toch(gsano[index12, 'mouse_homolog_gene'])
gs2 = intersect(gs2, colnames(p21$counts))
gl[[i]] = gs2
tmp = as.matrix(p21$counts[, gs2]) %>% rowMeans
alls = cbind(alls, score)
}
colnames(alls) = c("M/G1" , "G2" , "G2/M" , "G1/S" , "S phase" )
a=pheatmap::pheatmap(alls)
clu = a$tree_row
trr=cutree(clu,k=3)
cluster = paste('cluster',trr,sep='')
names(cluster) = names(trr)
annot=data.frame('cluster'=cluster,'group'=groups[names(cluster)],'Dpp4'=p21$counts[names(cluster),'Dpp4'],row.names=names(cluster))
x1 = alls[a$tree_row$order,a$tree_col$order]
pheatmap::pheatmap(x1[,2:5],annotation_row = annot,cluster_rows = F, cluster_cols = F,show_rownames = F,height=4,width=4.2, annotation_colors = anoc,border_color = "NA")
library(Seurat)
# A list of cell cycle markers
s.genes <- gl$`S phase`
g2m.genes <- gl$`G2/M`
# Create our Seurat object and complete the initalization steps
marrow <- CreateSeuratObject(counts = dcount)
Warning: Feature names cannot have underscores ('_'), replacing with dashes ('-')
marrow <- NormalizeData(marrow)
Performing log-normalization
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
marrow <- FindVariableFeatures(marrow, selection.method = "vst")
Calculating gene variances
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating feature variances of standardized and clipped values
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
marrow <- ScaleData(marrow, features = rownames(marrow))
Centering and scaling data matrix
|
| | 0%
|
|== | 2%
|
|==== | 4%
|
|====== | 6%
|
|======== | 7%
|
|========== | 9%
|
|============ | 11%
|
|============== | 13%
|
|================ | 15%
|
|================== | 17%
|
|==================== | 19%
|
|====================== | 20%
|
|======================== | 22%
|
|========================== | 24%
|
|=========================== | 26%
|
|============================= | 28%
|
|=============================== | 30%
|
|================================= | 31%
|
|=================================== | 33%
|
|===================================== | 35%
|
|======================================= | 37%
|
|========================================= | 39%
|
|=========================================== | 41%
|
|============================================= | 43%
|
|=============================================== | 44%
|
|================================================= | 46%
|
|=================================================== | 48%
|
|===================================================== | 50%
|
|======================================================= | 52%
|
|========================================================= | 54%
|
|=========================================================== | 56%
|
|============================================================= | 57%
|
|=============================================================== | 59%
|
|================================================================= | 61%
|
|=================================================================== | 63%
|
|===================================================================== | 65%
|
|======================================================================= | 67%
|
|========================================================================= | 69%
|
|=========================================================================== | 70%
|
|============================================================================= | 72%
|
|=============================================================================== | 74%
|
|================================================================================ | 76%
|
|================================================================================== | 78%
|
|==================================================================================== | 80%
|
|====================================================================================== | 81%
|
|======================================================================================== | 83%
|
|========================================================================================== | 85%
|
|============================================================================================ | 87%
|
|============================================================================================== | 89%
|
|================================================================================================ | 91%
|
|================================================================================================== | 93%
|
|==================================================================================================== | 94%
|
|====================================================================================================== | 96%
|
|======================================================================================================== | 98%
|
|==========================================================================================================| 100%
marrow <- RunPCA(marrow, features = VariableFeatures(marrow), ndims.print = 1:10, nfeatures.print = 10)
Warning in irlba(A = t(x = object), nv = npcs, ...) :
You're computing too large a percentage of total singular values, use a standard svd instead.
PC_ 1
Positive: Pira2, AC161409.1, Pirb, Lilra6, Lgals3, Ahnak, Lsp1, Fabp5, Ms4a6c, Olfm1
Negative: Gstm1, Ms4a3, Mpo, Tuba4a, Elane, Pglyrp1, Nedd4, Prtn3, Trem3, Pxylp1
PC_ 2
Positive: Ddx4, Il12a, Lgals1, Pclaf, Cnn3, Ssbp2, Tpm4, Hes6, Pcna, Rrm2
Negative: Ly6c2, Gsr, Prtn3, Igsf6, Rab32, Lcn2, Gda, Hp, Grn, Cldn15
PC_ 3
Positive: Ccna2, Ccnb1, Tubb4b, Aurkb, Kif22, H2afx, Plk1, Kif11, Top2a, Ube2c
Negative: Gstm1, Srgn, Nfkbia, RP23-221F6.1, Arap3, Rnd2, RP23-118E21.11, Mpo, RP23-405M23.1, Znfx1
PC_ 4
Positive: Rplp1, Rps27a, Rpl39, Rps26, RP23-123N23.5, Rpsa, Fos, mt-Rnr2, mt-Co2, mt-Rnr1
Negative: Rsrp1, RP23-23B17.2, Zyx, Abhd4, Rab37, Ppp1cb, Klhl6, Rasal3, Gcnt7, RP23-364M24.2
PC_ 5
Positive: Irf8, Rnf123, Rab19, Zfp36, Klf2, RP23-465M17.1, Arhgef6, Zfp68, Gnpda1, Sh2d3c
Negative: RP23-3F1.4, S100a4, Tlr13, Shpk, Zbed6, Tubb2a, P2ry6, Prdm9, Iqcc, Tmem106a
PC_ 6
Positive: Mcm5, Uhrf1, Mcm4, RP23-465M17.2, Pa2g4, Orc1, RP23-92B18.4, Slc35b4, Cdc14b, Calr
Negative: RP23-106P6.3, RP24-309H3.1, Chil3, Apip, Arl6ip1, RP23-405M23.1, Tbc1d12, Cd300d, H1f0, Mgst1
PC_ 7
Positive: Prr3, Braf, Ric1, Greb1, RP23-143O13.1, RP23-106I21.1, RP23-123N23.5, Zfp784, RP24-297N9.6, Tlr4
Negative: RP24-178J22.3, Pqlc2, Zfp119b, Lmln, RP23-456K21.1, Rhoq, RP23-456B18.2, Lace1, Adat1, RP23-221F6.1
PC_ 8
Positive: Vmp1, Ccdc91, Triqk, RP23-235P16.4, B3gnt3, Enox1, RP23-385F8.3, Syt14, Nell1, Cpn1
Negative: Skil, RP23-471B6.5, RP23-423B11.2, RP23-336M16.3, RP23-206J9.3, Klhl11, Il10ra, Zfp14, Crybg3, Igkv14-111
PC_ 9
Positive: Wnk4, Mllt4, Nfkbia, Zfp458, Cxcl2, Tmem94, AC123856.5, Plekhg5, Txnip, Apip
Negative: Dnajc16, Mib2, RP23-437N23.3, Lacc1, Phf13, Sgpl1, Rab40c, RP23-124N16.6, Ndrg1, Gtsf1
PC_ 10
Positive: Skil, Crybg3, RP23-471B6.5, RP23-358B23.3, Neb, Zfp595, RP23-32C18.6, Ccdc142os, Zfp182, Il6st
Negative: Zfp784, Cd33, Hpse, Hes6, Sidt2, Rad51c, Golga2, Slco5a1, Tbc1d12, Zfp942
marrow <- CellCycleScoring(marrow, s.features = s.genes, g2m.features = g2m.genes, set.ident = TRUE)
marrow <- RunUMAP(marrow, dims = 1:5)
Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
This message will be shown once per session
21:16:33 UMAP embedding parameters a = 0.9922 b = 1.112
21:16:33 Read 84 rows and found 5 numeric columns
21:16:33 Using Annoy for neighbor search, n_neighbors = 30
21:16:33 Building Annoy index with metric = cosine, n_trees = 50
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
21:16:33 Writing NN index file to temp file /tmp/RtmpQepCrE/file19d54709c48d
21:16:33 Searching Annoy index using 1 thread, search_k = 3000
21:16:33 Annoy recall = 100%
21:16:35 Commencing smooth kNN distance calibration using 1 thread
21:16:38 Initializing from normalized Laplacian + noise
21:16:39 Commencing optimization for 500 epochs, with 2618 positive edges
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
21:16:40 Optimization finished
DimPlot(marrow, reduction = "umap")
cluster = marrow@meta.data$seurat_clusters
names(cluster) = rownames(marrow@meta.data)
library(dplyr)
pbmc.markers <- FindAllMarkers(marrow, only.pos = TRUE, min.pct = 0.25, logfc.threshold = 0.25)
Calculating cluster 0
| | 0 % ~calculating
|+ | 1 % ~00s
|++ | 2 % ~00s
|++ | 3 % ~00s
|+++ | 4 % ~00s
|+++ | 5 % ~00s
|++++ | 7 % ~00s
|++++ | 8 % ~00s
|+++++ | 9 % ~00s
|+++++ | 10% ~00s
|++++++ | 11% ~00s
|+++++++ | 12% ~00s
|+++++++ | 13% ~00s
|++++++++ | 14% ~00s
|++++++++ | 15% ~00s
|+++++++++ | 16% ~00s
|+++++++++ | 18% ~00s
|++++++++++ | 19% ~00s
|++++++++++ | 20% ~00s
|+++++++++++ | 21% ~00s
|+++++++++++ | 22% ~00s
|++++++++++++ | 23% ~00s
|+++++++++++++ | 24% ~00s
|+++++++++++++ | 25% ~00s
|++++++++++++++ | 26% ~00s
|++++++++++++++ | 27% ~00s
|+++++++++++++++ | 29% ~00s
|+++++++++++++++ | 30% ~00s
|++++++++++++++++ | 31% ~00s
|++++++++++++++++ | 32% ~00s
|+++++++++++++++++ | 33% ~00s
|++++++++++++++++++ | 34% ~00s
|++++++++++++++++++ | 35% ~00s
|+++++++++++++++++++ | 36% ~00s
|+++++++++++++++++++ | 37% ~00s
|++++++++++++++++++++ | 38% ~00s
|++++++++++++++++++++ | 40% ~00s
|+++++++++++++++++++++ | 41% ~00s
|+++++++++++++++++++++ | 42% ~00s
|++++++++++++++++++++++ | 43% ~00s
|++++++++++++++++++++++ | 44% ~00s
|+++++++++++++++++++++++ | 45% ~00s
|++++++++++++++++++++++++ | 46% ~00s
|++++++++++++++++++++++++ | 47% ~00s
|+++++++++++++++++++++++++ | 48% ~00s
|+++++++++++++++++++++++++ | 49% ~00s
|++++++++++++++++++++++++++ | 51% ~00s
|++++++++++++++++++++++++++ | 52% ~00s
|+++++++++++++++++++++++++++ | 53% ~00s
|+++++++++++++++++++++++++++ | 54% ~00s
|++++++++++++++++++++++++++++ | 55% ~00s
|+++++++++++++++++++++++++++++ | 56% ~00s
|+++++++++++++++++++++++++++++ | 57% ~00s
|++++++++++++++++++++++++++++++ | 58% ~00s
|++++++++++++++++++++++++++++++ | 59% ~00s
|+++++++++++++++++++++++++++++++ | 60% ~00s
|+++++++++++++++++++++++++++++++ | 62% ~00s
|++++++++++++++++++++++++++++++++ | 63% ~00s
|++++++++++++++++++++++++++++++++ | 64% ~00s
|+++++++++++++++++++++++++++++++++ | 65% ~00s
|+++++++++++++++++++++++++++++++++ | 66% ~00s
|++++++++++++++++++++++++++++++++++ | 67% ~00s
|+++++++++++++++++++++++++++++++++++ | 68% ~00s
|+++++++++++++++++++++++++++++++++++ | 69% ~00s
|++++++++++++++++++++++++++++++++++++ | 70% ~00s
|++++++++++++++++++++++++++++++++++++ | 71% ~00s
|+++++++++++++++++++++++++++++++++++++ | 73% ~00s
|+++++++++++++++++++++++++++++++++++++ | 74% ~00s
|++++++++++++++++++++++++++++++++++++++ | 75% ~00s
|++++++++++++++++++++++++++++++++++++++ | 76% ~00s
|+++++++++++++++++++++++++++++++++++++++ | 77% ~00s
|++++++++++++++++++++++++++++++++++++++++ | 78% ~00s
|++++++++++++++++++++++++++++++++++++++++ | 79% ~00s
|+++++++++++++++++++++++++++++++++++++++++ | 80% ~00s
|+++++++++++++++++++++++++++++++++++++++++ | 81% ~00s
|++++++++++++++++++++++++++++++++++++++++++ | 82% ~00s
|++++++++++++++++++++++++++++++++++++++++++ | 84% ~00s
|+++++++++++++++++++++++++++++++++++++++++++ | 85% ~00s
|+++++++++++++++++++++++++++++++++++++++++++ | 86% ~00s
|++++++++++++++++++++++++++++++++++++++++++++ | 87% ~00s
|++++++++++++++++++++++++++++++++++++++++++++ | 88% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++ | 89% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++ | 90% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++ | 91% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++ | 92% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++ | 93% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++ | 95% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++ | 96% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++++ | 97% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++++ | 98% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++++| 99% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=00s
Calculating cluster 1
| | 0 % ~calculating
|+ | 1 % ~00s
|++ | 2 % ~00s
|++ | 3 % ~00s
|+++ | 4 % ~00s
|+++ | 5 % ~00s
|++++ | 6 % ~00s
|++++ | 7 % ~00s
|+++++ | 8 % ~00s
|+++++ | 9 % ~00s
|++++++ | 10% ~00s
|++++++ | 11% ~00s
|+++++++ | 12% ~00s
|+++++++ | 13% ~00s
|++++++++ | 14% ~00s
|++++++++ | 15% ~00s
|+++++++++ | 16% ~00s
|+++++++++ | 17% ~00s
|++++++++++ | 18% ~00s
|++++++++++ | 19% ~00s
|+++++++++++ | 20% ~00s
|+++++++++++ | 21% ~00s
|++++++++++++ | 22% ~00s
|++++++++++++ | 23% ~00s
|+++++++++++++ | 24% ~00s
|+++++++++++++ | 26% ~00s
|++++++++++++++ | 27% ~00s
|++++++++++++++ | 28% ~00s
|+++++++++++++++ | 29% ~00s
|+++++++++++++++ | 30% ~00s
|++++++++++++++++ | 31% ~00s
|++++++++++++++++ | 32% ~00s
|+++++++++++++++++ | 33% ~00s
|+++++++++++++++++ | 34% ~00s
|++++++++++++++++++ | 35% ~00s
|++++++++++++++++++ | 36% ~00s
|+++++++++++++++++++ | 37% ~00s
|+++++++++++++++++++ | 38% ~00s
|++++++++++++++++++++ | 39% ~00s
|++++++++++++++++++++ | 40% ~00s
|+++++++++++++++++++++ | 41% ~00s
|+++++++++++++++++++++ | 42% ~00s
|++++++++++++++++++++++ | 43% ~00s
|++++++++++++++++++++++ | 44% ~00s
|+++++++++++++++++++++++ | 45% ~00s
|+++++++++++++++++++++++ | 46% ~00s
|++++++++++++++++++++++++ | 47% ~00s
|++++++++++++++++++++++++ | 48% ~00s
|+++++++++++++++++++++++++ | 49% ~00s
|+++++++++++++++++++++++++ | 50% ~00s
|++++++++++++++++++++++++++ | 51% ~00s
|+++++++++++++++++++++++++++ | 52% ~00s
|+++++++++++++++++++++++++++ | 53% ~00s
|++++++++++++++++++++++++++++ | 54% ~00s
|++++++++++++++++++++++++++++ | 55% ~00s
|+++++++++++++++++++++++++++++ | 56% ~00s
|+++++++++++++++++++++++++++++ | 57% ~00s
|++++++++++++++++++++++++++++++ | 58% ~00s
|++++++++++++++++++++++++++++++ | 59% ~00s
|+++++++++++++++++++++++++++++++ | 60% ~00s
|+++++++++++++++++++++++++++++++ | 61% ~00s
|++++++++++++++++++++++++++++++++ | 62% ~00s
|++++++++++++++++++++++++++++++++ | 63% ~00s
|+++++++++++++++++++++++++++++++++ | 64% ~00s
|+++++++++++++++++++++++++++++++++ | 65% ~00s
|++++++++++++++++++++++++++++++++++ | 66% ~00s
|++++++++++++++++++++++++++++++++++ | 67% ~00s
|+++++++++++++++++++++++++++++++++++ | 68% ~00s
|+++++++++++++++++++++++++++++++++++ | 69% ~00s
|++++++++++++++++++++++++++++++++++++ | 70% ~00s
|++++++++++++++++++++++++++++++++++++ | 71% ~00s
|+++++++++++++++++++++++++++++++++++++ | 72% ~00s
|+++++++++++++++++++++++++++++++++++++ | 73% ~00s
|++++++++++++++++++++++++++++++++++++++ | 74% ~00s
|++++++++++++++++++++++++++++++++++++++ | 76% ~00s
|+++++++++++++++++++++++++++++++++++++++ | 77% ~00s
|+++++++++++++++++++++++++++++++++++++++ | 78% ~00s
|++++++++++++++++++++++++++++++++++++++++ | 79% ~00s
|++++++++++++++++++++++++++++++++++++++++ | 80% ~00s
|+++++++++++++++++++++++++++++++++++++++++ | 81% ~00s
|+++++++++++++++++++++++++++++++++++++++++ | 82% ~00s
|++++++++++++++++++++++++++++++++++++++++++ | 83% ~00s
|++++++++++++++++++++++++++++++++++++++++++ | 84% ~00s
|+++++++++++++++++++++++++++++++++++++++++++ | 85% ~00s
|+++++++++++++++++++++++++++++++++++++++++++ | 86% ~00s
|++++++++++++++++++++++++++++++++++++++++++++ | 87% ~00s
|++++++++++++++++++++++++++++++++++++++++++++ | 88% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++ | 89% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++ | 90% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++ | 91% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++ | 92% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++ | 93% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++ | 94% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++ | 95% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++ | 96% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++++ | 97% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++++ | 98% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++++| 99% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=00s
Calculating cluster 2
| | 0 % ~calculating
|+ | 1 % ~00s
|++ | 2 % ~00s
|++ | 3 % ~00s
|+++ | 4 % ~00s
|+++ | 5 % ~00s
|++++ | 7 % ~00s
|++++ | 8 % ~00s
|+++++ | 9 % ~00s
|+++++ | 10% ~00s
|++++++ | 11% ~00s
|++++++ | 12% ~00s
|+++++++ | 13% ~00s
|++++++++ | 14% ~00s
|++++++++ | 15% ~00s
|+++++++++ | 16% ~00s
|+++++++++ | 17% ~00s
|++++++++++ | 18% ~00s
|++++++++++ | 20% ~00s
|+++++++++++ | 21% ~00s
|+++++++++++ | 22% ~00s
|++++++++++++ | 23% ~00s
|++++++++++++ | 24% ~00s
|+++++++++++++ | 25% ~00s
|++++++++++++++ | 26% ~00s
|++++++++++++++ | 27% ~00s
|+++++++++++++++ | 28% ~00s
|+++++++++++++++ | 29% ~00s
|++++++++++++++++ | 30% ~00s
|++++++++++++++++ | 32% ~00s
|+++++++++++++++++ | 33% ~00s
|+++++++++++++++++ | 34% ~00s
|++++++++++++++++++ | 35% ~00s
|++++++++++++++++++ | 36% ~00s
|+++++++++++++++++++ | 37% ~00s
|++++++++++++++++++++ | 38% ~00s
|++++++++++++++++++++ | 39% ~00s
|+++++++++++++++++++++ | 40% ~00s
|+++++++++++++++++++++ | 41% ~00s
|++++++++++++++++++++++ | 42% ~00s
|++++++++++++++++++++++ | 43% ~00s
|+++++++++++++++++++++++ | 45% ~00s
|+++++++++++++++++++++++ | 46% ~00s
|++++++++++++++++++++++++ | 47% ~00s
|++++++++++++++++++++++++ | 48% ~00s
|+++++++++++++++++++++++++ | 49% ~00s
|+++++++++++++++++++++++++ | 50% ~00s
|++++++++++++++++++++++++++ | 51% ~00s
|+++++++++++++++++++++++++++ | 52% ~00s
|+++++++++++++++++++++++++++ | 53% ~00s
|++++++++++++++++++++++++++++ | 54% ~00s
|++++++++++++++++++++++++++++ | 55% ~00s
|+++++++++++++++++++++++++++++ | 57% ~00s
|+++++++++++++++++++++++++++++ | 58% ~00s
|++++++++++++++++++++++++++++++ | 59% ~00s
|++++++++++++++++++++++++++++++ | 60% ~00s
|+++++++++++++++++++++++++++++++ | 61% ~00s
|+++++++++++++++++++++++++++++++ | 62% ~00s
|++++++++++++++++++++++++++++++++ | 63% ~00s
|+++++++++++++++++++++++++++++++++ | 64% ~00s
|+++++++++++++++++++++++++++++++++ | 65% ~00s
|++++++++++++++++++++++++++++++++++ | 66% ~00s
|++++++++++++++++++++++++++++++++++ | 67% ~00s
|+++++++++++++++++++++++++++++++++++ | 68% ~00s
|+++++++++++++++++++++++++++++++++++ | 70% ~00s
|++++++++++++++++++++++++++++++++++++ | 71% ~00s
|++++++++++++++++++++++++++++++++++++ | 72% ~00s
|+++++++++++++++++++++++++++++++++++++ | 73% ~00s
|+++++++++++++++++++++++++++++++++++++ | 74% ~00s
|++++++++++++++++++++++++++++++++++++++ | 75% ~00s
|+++++++++++++++++++++++++++++++++++++++ | 76% ~00s
|+++++++++++++++++++++++++++++++++++++++ | 77% ~00s
|++++++++++++++++++++++++++++++++++++++++ | 78% ~00s
|++++++++++++++++++++++++++++++++++++++++ | 79% ~00s
|+++++++++++++++++++++++++++++++++++++++++ | 80% ~00s
|+++++++++++++++++++++++++++++++++++++++++ | 82% ~00s
|++++++++++++++++++++++++++++++++++++++++++ | 83% ~00s
|++++++++++++++++++++++++++++++++++++++++++ | 84% ~00s
|+++++++++++++++++++++++++++++++++++++++++++ | 85% ~00s
|+++++++++++++++++++++++++++++++++++++++++++ | 86% ~00s
|++++++++++++++++++++++++++++++++++++++++++++ | 87% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++ | 88% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++ | 89% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++ | 90% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++ | 91% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++ | 92% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++ | 93% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++ | 95% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++ | 96% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++++ | 97% ~00s
|+++++++++++++++++++++++++++++++++++++++++++++++++ | 98% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++++| 99% ~00s
|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=00s
pbmc.markers %>% group_by(cluster) %>% top_n(n = 2, wt = avg_log2FC)
pbmc.markers =pbmc.markers[pbmc.markers$p_val_adj<0.05,]
top10 <- pbmc.markers %>% group_by(cluster) %>% top_n(n = 10, wt = avg_log2FC)
#write.csv(top10,'top20.csv',quote=F)
top upregulated genes
p = DoHeatmap(pbmc, features = top10$gene,group.colors = c('green4','royalblue4','goldenrod1'))
#ggsave('cluster.seurate.pdf',p,height=8,width=5)
p
exp = pbmc@assays$RNA@data
a1=vlnplot('Dpp4',cluster,exp,colp=NULL,pt.size=0.1,height=0.6)
Scale for 'y' is already present. Adding another scale for 'y', which will replace the existing scale.
a1+ ggpubr::stat_compare_means(label = "p.signif", label.x = 1.5)
source('../lib.r')
treat = 'P'
ref='NP'
group = groups
group=group[group %in% c(treat,ref)]
group[group!=ref]='treat'
group[group==ref]='control'
res1=runDE(dcount,group,ref='control')
group
control treat
33 34
converting counts to integer mode
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
-- replacing outliers and refitting for 4021 genes
-- DESeq argument 'minReplicatesForReplace' = 7
-- original counts are preserved in counts(dds)
estimating dispersions
fitting model and testing
res1[1:5,]
DotPlot.DE = function(de1,num=10,ylab='Treat vs control',gn1=NULL,gn2=NULL,orderpvalue=NULL,fc=2,pval=NULL){
library(ggplot2)
library(ggrepel)
tmp=de1
tmp=tmp[!is.na(tmp$pvalue),]
tmp$score=-log10(tmp$padj)
if (!is.null(pval)){
tmp$score=-log10(tmp$pvalue)
}
gg=data.frame('log2FoldChange'=tmp$log2FoldChange,'score'=tmp$score,'Z'=tmp$Z,'name'=rownames(tmp),
'padj'=tmp$padj,'pvalue'=tmp$pvalue,'score2'=tmp$score*sign(tmp$log2FoldChange))
p=ggplot(gg,aes(y=score,x=log2FoldChange))
p=p+geom_point(shape=".",alpha=1/1,size = 1,color='grey')+theme_classic()
p=p+labs(y=ylab)
gg$score2=sign(gg$log2FoldChange)*gg$score
rownames(gg)=gg$name
if (is.null(orderpvalue)){
t1=gg[order(gg$log2FoldChange),]
t2=gg[order(gg$log2FoldChange,decreasing=TRUE),]
}else{
t1=gg[order(gg$score2),]
t2=gg[order(gg$score2,decreasing=TRUE),]
}
t1=t1[t1$log2FoldChange<(-1),]
n1=rownames(t1)[1:num]
t2=t2[t2$log2FoldChange>1,]
n2=rownames(t2)[1:num]
if (!is.null(gn1)){
n1=c(n1,gn1)
}
if (!is.null(gn2)){
n1=c(n2,gn2)
}
p=p+geom_text_repel(
data = gg[gg$name %in% n1,],
aes(label = name),
size = 3,
box.padding = unit(0.35, "lines"),
point.padding = unit(0.3, "lines"),
color='blue'
)
p=p+geom_text_repel(
data = gg[gg$name %in% n2,],
aes(label = name),
size = 3,
box.padding = unit(0.35, "lines"),
point.padding = unit(0.3, "lines"),
color='red1'
)
p=p+ geom_point(data = subset(gg,(log2FoldChange< (-fc) & pvalue < 0.01)),shape=".",alpha=1/1,size = 2,color='blue')
p=p+ geom_point(data = subset(gg,(log2FoldChange>fc & pvalue < 0.01)),shape=".",alpha=1/1,size = 2,color='red1')
return(p)
}
DotPlot.DE(res1,pval=TRUE,ylab='-10log padj')
Warning: ggrepel: 8 unlabeled data points (too many overlaps). Consider increasing max.overlaps
ref='NP'
treat = 'P-PNP'
group = groups
group[group %in% c('P','PNP')]='P-PNP'
table(group)
group
NP P-PNP
33 51
group=group[group %in% c(treat,ref)]
group[group!=ref]='treat'
group[group==ref]='control'
res1=runDE(dcount,group,ref='control')
group
control treat
33 51
converting counts to integer mode
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
-- replacing outliers and refitting for 4354 genes
-- DESeq argument 'minReplicatesForReplace' = 7
-- original counts are preserved in counts(dds)
estimating dispersions
fitting model and testing
res1[1:5,]