From b3ec860e0b7ee8dbde319e29774e855174a2acca Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 12 Mar 2021 09:44:47 -0500
Subject: [PATCH 01/74] Update README.md
---
README.md | 1 +
1 file changed, 1 insertion(+)
diff --git a/README.md b/README.md
index ad4361e..15f93a1 100644
--- a/README.md
+++ b/README.md
@@ -15,6 +15,7 @@ iCellR is an interactive R package to work with high-throughput single cell sequ
- Tutorial: [example 1 code](https://genome.med.nyu.edu/results/external/iCellR/example1/code.txt) and [results](https://genome.med.nyu.edu/results/external/iCellR/example1/) (based on KNetL map )
- Tutorial: [example 2 code](https://genome.med.nyu.edu/results/external/iCellR/example2/code.txt) and [results](https://genome.med.nyu.edu/results/external/iCellR/example2/) (based on CPCA batch alignment and KNetL map )
- Link to a video tutorial for CITE-Seq and scRNA-Seq analysis: [Video](https://vimeo.com/337822487)
+- All you need to know about KNetL map: [Video](https://youtu.be/tkoPTVciQm0)
- Link to manual [Manual](https://cran.r-project.org/web/packages/iCellR/iCellR.pdf) and Comprehensive R Archive Network [(CRAN)](https://cran.r-project.org/web/packages/iCellR/index.html).
iCellR Viewer (web GUI app): https://compbio.nyumc.org/icellr/
From 09c15f853dd56a340a0fad7f3070209edf7a0adb Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 12 Mar 2021 09:49:12 -0500
Subject: [PATCH 02/74] Update README.md
---
README.md | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/README.md b/README.md
index 15f93a1..7e66b8b 100644
--- a/README.md
+++ b/README.md
@@ -422,9 +422,9 @@ my.obj <- run.pc.tsne(my.obj, dims = 1:10)
# UMAP
my.obj <- run.umap(my.obj, dims = 1:10)
-# KNetL (for lager than 5000 cell use a k of about 400)
+# KNetL (for lager than 5000 cell use a zoom of about 400)
# Because knetl has a very high resolution it's best to use a dim of 20 (this usually works best for most data)
-my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap")
+my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap") # (Important note!) don't forget to set the zoom in the right range
########################### IMPORTANT NOTE ########################################
#### Because KNetl has a very high resolution it's best to use a dim of 20 (this usually works best for most data)
From 6582acf6670488f6b23f860331177c2d02684bcc Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 12 Mar 2021 10:15:32 -0500
Subject: [PATCH 03/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 7e66b8b..3b7910d 100644
--- a/README.md
+++ b/README.md
@@ -411,7 +411,7 @@ opt.pcs.plot(my.obj)
-- Perform other dimensionality reductiond (tSNE, UMAP, KNetL, PHATE, diffusion map)
+- Perform other dimensionality reductions (tSNE, UMAP, KNetL, PHATE, destiny, diffusion maps)
We recommend tSNE, UMAP and KNetL. KNetL is fundamentally more powerful.
From 93c5a0a1bf85a1ab26cb907ff9a2879f09b203a4 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 12 Mar 2021 10:20:29 -0500
Subject: [PATCH 04/74] Update README.md
---
README.md | 4 ++++
1 file changed, 4 insertions(+)
diff --git a/README.md b/README.md
index 3b7910d..3ceb716 100644
--- a/README.md
+++ b/README.md
@@ -510,6 +510,10 @@ This is one of the harder parts of the analysis and sometimes you need to adjust
my.obj <- iclust(my.obj, sensitivity = 150, data.type = "knetl")
+# clustering based on PCA
+
+# my.obj <- iclust(my.obj, sensitivity = 150, data.type = "pca", dims=1:10)
+
# play with k to get the clusters right. Usually 150 is good.
###### more examples
From e1442280e52e3c3ca52547de017b9cd453c9232d Mon Sep 17 00:00:00 2001
From: Khodadadi-Jamayran
Date: Mon, 29 Mar 2021 19:03:02 -0400
Subject: [PATCH 05/74] ST
---
DESCRIPTION | 4 ++--
R/F0021.R | 19 +++++++++++++------
R/F0028.R | 3 +++
R/F0070.R | 26 +++++++++++++++-----------
man/clust.avg.exp.Rd | 12 +++++++++++-
man/spatial.plot.Rd | 2 +-
6 files changed, 45 insertions(+), 21 deletions(-)
diff --git a/DESCRIPTION b/DESCRIPTION
index 0f6f05d..99405d2 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
Package: iCellR
Type: Package
Title: Analyzing High-Throughput Single Cell Sequencing Data
-Version: 1.6.1
+Version: 1.6.2
Authors@R: c(
person(given = 'Alireza', family = 'Khodadadi-Jamayran',role = c('aut','cre'), email = 'alireza.khodadadi.j@gmail.com', comment = c(ORCID = '0000-0003-2495-7504')),
person(given = 'Joseph', family = 'Pucella',role = c('aut','ctb'), comment = c(ORCID = '0000-0003-0875-8046')),
@@ -27,7 +27,7 @@ RoxygenNote: 7.1.1
BugReports: https://github.com/rezakj/iCellR/issues
URL: https://github.com/rezakj/iCellR
NeedsCompilation: yes
-Packaged: 2021-03-03 22:59:24 UTC; khodaa01
+Packaged: 2021-03-29 23:01:15 UTC; khodaa01
Author: Alireza Khodadadi-Jamayran [aut, cre]
(),
Joseph Pucella [aut, ctb] (),
diff --git a/R/F0021.R b/R/F0021.R
index b090d1e..cbb1195 100644
--- a/R/F0021.R
+++ b/R/F0021.R
@@ -4,12 +4,16 @@
#' @param x An object of class iCellR.
#' @param data.type Choose from "main" and "imputed", default = "main"
#' @param conds.to.avg Choose the conditions you want to average, default = NULL (all conditions).
+#' @param rounding.digits integer indicating the number of decimal places (round) or significant digits (signif) to be used.
+#' @param round.num Rounding of Numbers, default = FALSE.
#' @return An object of class iCellR.
#' @import progress
#' @export
clust.avg.exp <- function (x = NULL,
data.type = "main",
- conds.to.avg = NULL) {
+ conds.to.avg = NULL,
+ rounding.digits = 4,
+ round.num = FALSE) {
if ("iCellR" != class(x)[1]) {
stop("x should be an object of class iCellR")
}
@@ -55,16 +59,17 @@ clust.avg.exp <- function (x = NULL,
colnames(DATA) <- NameCol
DATA <- cbind(gene = rownames(DATA), DATA)
rownames(DATA) <- NULL
-# eval(call("<-", as.name(NameCol), DATA))
- datalist[[i]] <- DATA
+ eval(call("<-", as.name(NameCol), DATA))
+# datalist[[i]] <- DATA
}
# multmerge = function(mypath){
# filenames=list.files(pattern="meanExp")
# datalist = lapply(filenames, function(x){read.table(file=x,header=T)})
# Reduce(function(x,y) {merge(x,y)}, datalist)
# }
-# filenames <- ls(pattern="cluster_")
-# datalist <- mget(filenames)
+ filenames <- ls(pattern="cluster_")
+ filenames <- filenames[order(nchar(filenames))]
+ datalist <- mget(filenames)
MeanExpForClusters <- Reduce(function(x,y) {merge(x,y)}, datalist)
#
# MeanExpForClusters <- multmerge()
@@ -79,7 +84,9 @@ clust.avg.exp <- function (x = NULL,
row.names(data) <- data$gene
data <- data[,-1]
# data = as.data.frame(sapply(data, as.numeric))
- data <- round(data,digits=4)
+ if (round.num == TRUE) {
+ data <- round(data, digits = rounding.digits)
+ }
data <- cbind(gene=MeanExpForClusters$gene,data)
# MeanExpForClusters <- MeanExpForClusters[order(nchar(colnames(MeanExpForClusters)),colnames(MeanExpForClusters))]
attributes(x)$clust.avg <- data
diff --git a/R/F0028.R b/R/F0028.R
index 3012058..784f56c 100644
--- a/R/F0028.R
+++ b/R/F0028.R
@@ -38,6 +38,9 @@ findMarkers <- function (x = NULL,
# get avrages
x <- clust.avg.exp(x, data.type = data.type)
DATA <- x@best.clust
+ if(!is.numeric(DATA$clusters)){
+ stop("Cluster names have to be numeric")
+ }
############## set wich clusters you want as condition 1 and 2
MyClusts <- as.numeric(unique(DATA$clusters))
MyClusts <- sort(MyClusts)
diff --git a/R/F0070.R b/R/F0070.R
index c64c351..8ef43b4 100644
--- a/R/F0070.R
+++ b/R/F0070.R
@@ -31,7 +31,7 @@
#' @importFrom ggplot2 ggplot theme_classic geom_segment geom_violin guide_colorbar guide_legend guides scale_color_discrete scale_colour_gradient scale_fill_gradient2 scale_x_continuous scale_y_continuous scale_y_discrete stat_summary coord_polar element_rect element_text element_blank facet_wrap scale_color_manual geom_hline geom_jitter geom_vline ylab xlab ggtitle theme_bw aes theme geom_bar geom_point geom_boxplot geom_errorbar position_dodge geom_tile geom_density geom_line
#' @export
spatial.plot <- function (x = NULL,
- cell.size = 0.5,
+ cell.size = 1,
cell.colors = c("gray","red"),
back.col = "black",
col.by = "clusters",
@@ -164,6 +164,8 @@ if (col.by == "cc") {
if (!is.null(conds.to.plot)) {
data <- subset(data, data$MYconds %in% conds.to.plot)
}
+################# Fix V5
+# data$V5 <- data$V5 + ((min(data$V5)) * -1)
############################# plot
if (!is.numeric(data$clusters)) {
if(length(MYConds) == 1) {
@@ -185,8 +187,10 @@ if (col.by == "cc") {
############# Annotation
if (anno.clust == TRUE) {
cords <- aggregate(data[, 6:7], list(data$clusters), mean)
- MYX=(cords$V5) * -1
- MYY=(cords$V6) * -1
+# MYX=(cords$V5) * -1
+# MYY=(cords$V6) * -1
+ MYX=(cords$V6)
+ MYY=(cords$V5)
MYZ=cords$Group.1
myPLOT <- myPLOT + annotate("text",
x = MYX,
@@ -204,9 +208,9 @@ if (col.by == "cc") {
guides(colour = guide_legend(override.aes = list(size=5))) +
scale_color_discrete(name="") +
ggtitle(MyTitle) +
- theme(panel.background = element_rect(fill = "black", colour = "black"),
+ theme(panel.background = element_rect(fill = back.col, colour = back.col),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
- legend.key = element_rect(fill = "black")) +
+ legend.key = element_rect(fill = back.col)) +
theme(axis.title.x=element_blank(),
axis.text.x=element_blank(),
axis.ticks.x=element_blank()) +
@@ -216,8 +220,8 @@ if (col.by == "cc") {
############# Annotation
if (anno.clust == TRUE) {
cords <- aggregate(data[, 6:7], list(data$clusters), mean)
- MYX=(cords$V5) * -1
- MYY=(cords$V6) * -1
+ MYX=(cords$V6)
+ MYY=(cords$V5)
MYZ=cords$Group.1
myPLOT <- myPLOT + annotate("text",
x = MYX,
@@ -236,9 +240,9 @@ if (col.by == "cc") {
guides(colour = guide_legend(override.aes = list(size=5))) +
scale_color_discrete(name="") +
ggtitle(MyTitle) +
- theme(panel.background = element_rect(fill = "black", colour = "black"),
+ theme(panel.background = element_rect(fill = back.col, colour = back.col),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
- legend.key = element_rect(fill = "black")) +
+ legend.key = element_rect(fill = back.col)) +
theme(axis.title.x=element_blank(),
axis.text.x=element_blank(),
axis.ticks.x=element_blank()) +
@@ -274,9 +278,9 @@ if (col.by == "cc") {
geom_point(size = cell.size, alpha = cell.transparency) +
scale_colour_gradient(low = cell.colors[1], high = cell.colors[2], name="") +
ggtitle(MyTitle) +
- theme(panel.background = element_rect(fill = "black", colour = "black"),
+ theme(panel.background = element_rect(fill = back.col, colour = back.col),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
- legend.key = element_rect(fill = "black")) +
+ legend.key = element_rect(fill = back.col)) +
theme(axis.title.x=element_blank(),
axis.text.x=element_blank(),
axis.ticks.x=element_blank()) +
diff --git a/man/clust.avg.exp.Rd b/man/clust.avg.exp.Rd
index 13ab719..97c58c5 100644
--- a/man/clust.avg.exp.Rd
+++ b/man/clust.avg.exp.Rd
@@ -4,7 +4,13 @@
\alias{clust.avg.exp}
\title{Create a data frame of mean expression of genes per cluster}
\usage{
-clust.avg.exp(x = NULL, data.type = "main", conds.to.avg = NULL)
+clust.avg.exp(
+ x = NULL,
+ data.type = "main",
+ conds.to.avg = NULL,
+ rounding.digits = 4,
+ round.num = FALSE
+)
}
\arguments{
\item{x}{An object of class iCellR.}
@@ -12,6 +18,10 @@ clust.avg.exp(x = NULL, data.type = "main", conds.to.avg = NULL)
\item{data.type}{Choose from "main" and "imputed", default = "main"}
\item{conds.to.avg}{Choose the conditions you want to average, default = NULL (all conditions).}
+
+\item{rounding.digits}{integer indicating the number of decimal places (round) or significant digits (signif) to be used.}
+
+\item{round.num}{Rounding of Numbers, default = FALSE.}
}
\value{
An object of class iCellR.
diff --git a/man/spatial.plot.Rd b/man/spatial.plot.Rd
index 92dbf59..d10df7d 100644
--- a/man/spatial.plot.Rd
+++ b/man/spatial.plot.Rd
@@ -6,7 +6,7 @@
\usage{
spatial.plot(
x = NULL,
- cell.size = 0.5,
+ cell.size = 1,
cell.colors = c("gray", "red"),
back.col = "black",
col.by = "clusters",
From 89c4f28022cdbe43a24eb49e3f907cd319d41ea4 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 29 Mar 2021 19:05:44 -0400
Subject: [PATCH 06/74] Update README.md
---
README.md | 11 ++++++++---
1 file changed, 8 insertions(+), 3 deletions(-)
diff --git a/README.md b/README.md
index ad4361e..3ceb716 100644
--- a/README.md
+++ b/README.md
@@ -15,6 +15,7 @@ iCellR is an interactive R package to work with high-throughput single cell sequ
- Tutorial: [example 1 code](https://genome.med.nyu.edu/results/external/iCellR/example1/code.txt) and [results](https://genome.med.nyu.edu/results/external/iCellR/example1/) (based on KNetL map )
- Tutorial: [example 2 code](https://genome.med.nyu.edu/results/external/iCellR/example2/code.txt) and [results](https://genome.med.nyu.edu/results/external/iCellR/example2/) (based on CPCA batch alignment and KNetL map )
- Link to a video tutorial for CITE-Seq and scRNA-Seq analysis: [Video](https://vimeo.com/337822487)
+- All you need to know about KNetL map: [Video](https://youtu.be/tkoPTVciQm0)
- Link to manual [Manual](https://cran.r-project.org/web/packages/iCellR/iCellR.pdf) and Comprehensive R Archive Network [(CRAN)](https://cran.r-project.org/web/packages/iCellR/index.html).
iCellR Viewer (web GUI app): https://compbio.nyumc.org/icellr/
@@ -410,7 +411,7 @@ opt.pcs.plot(my.obj)
-- Perform other dimensionality reductiond (tSNE, UMAP, KNetL, PHATE, diffusion map)
+- Perform other dimensionality reductions (tSNE, UMAP, KNetL, PHATE, destiny, diffusion maps)
We recommend tSNE, UMAP and KNetL. KNetL is fundamentally more powerful.
@@ -421,9 +422,9 @@ my.obj <- run.pc.tsne(my.obj, dims = 1:10)
# UMAP
my.obj <- run.umap(my.obj, dims = 1:10)
-# KNetL (for lager than 5000 cell use a k of about 400)
+# KNetL (for lager than 5000 cell use a zoom of about 400)
# Because knetl has a very high resolution it's best to use a dim of 20 (this usually works best for most data)
-my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap")
+my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap") # (Important note!) don't forget to set the zoom in the right range
########################### IMPORTANT NOTE ########################################
#### Because KNetl has a very high resolution it's best to use a dim of 20 (this usually works best for most data)
@@ -509,6 +510,10 @@ This is one of the harder parts of the analysis and sometimes you need to adjust
my.obj <- iclust(my.obj, sensitivity = 150, data.type = "knetl")
+# clustering based on PCA
+
+# my.obj <- iclust(my.obj, sensitivity = 150, data.type = "pca", dims=1:10)
+
# play with k to get the clusters right. Usually 150 is good.
###### more examples
From 430e783115297fb52c64efb53eef098f9dff915d Mon Sep 17 00:00:00 2001
From: Khodadadi-Jamayran
Date: Thu, 1 Apr 2021 13:48:45 -0400
Subject: [PATCH 07/74] baseMeanFix
---
DESCRIPTION | 2 +-
R/F0028.R | 4 ++--
2 files changed, 3 insertions(+), 3 deletions(-)
diff --git a/DESCRIPTION b/DESCRIPTION
index 99405d2..cbef2b6 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -27,7 +27,7 @@ RoxygenNote: 7.1.1
BugReports: https://github.com/rezakj/iCellR/issues
URL: https://github.com/rezakj/iCellR
NeedsCompilation: yes
-Packaged: 2021-03-29 23:01:15 UTC; khodaa01
+Packaged: 2021-04-01 17:47:14 UTC; khodaa01
Author: Alireza Khodadadi-Jamayran [aut, cre]
(),
Joseph Pucella [aut, ctb] (),
diff --git a/R/F0028.R b/R/F0028.R
index 784f56c..08237f6 100644
--- a/R/F0028.R
+++ b/R/F0028.R
@@ -172,8 +172,8 @@ findMarkers <- function (x = NULL,
dfm <- dfm[order(myClustOrd, decreasing = FALSE),]
df <- as.data.frame(dfm)
#######
- df$clusters <- as.numeric(df$clusters)
- df$baseMean <- as.numeric(df$baseMean)
+ df$clusters <- as.numeric(as.character(df$clusters))
+ df$baseMean <- as.numeric(as.character(df$baseMean))
######
message("All done!")
return(df)
From 13befd5a4c2a6e26f2182312eb30d0b943636259 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 10:57:21 -0400
Subject: [PATCH 08/74] Update README.md
---
README.md | 142 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 142 insertions(+)
diff --git a/README.md b/README.md
index 3ceb716..9dcd47e 100644
--- a/README.md
+++ b/README.md
@@ -2739,3 +2739,145 @@ dev.off()
+# Single cell ATAC sequencing with scRNA-Seq (scATAC-Seq)
+
+```r
+library("iCellR")
+my.data <- load10x("filtered_gene_bc_matrices/")
+
+# see the row names
+row.names(my.data)
+
+# get peak names
+ATAC <- grep("^chr",row.names(my.data),value=T)
+
+# get scATAC data
+MyATAC <- subset(my.data, row.names(my.data) %in% ATAC)
+head(MyATAC)[1:3]
+# AAACAGCCAAGTGAAC.1 AAACAGCCACTGACCG.1 AAACAGCCATGATTGT.1
+#chr1.181218.181695 0 0 1
+#chr1.191296.191699 0 0 0
+#chr1.629770.630129 0 0 0
+#chr1.633806.634251 0 0 0
+#chr1.778422.779040 0 0 0
+#chr1.827306.827702 0 0 0
+
+dim(MyATAC)
+# [1] 21923 6326
+
+# get RNA data
+MyRNAs <- subset(my.data, !row.names(my.data) %in% ATAC)
+head(MyRNAs)[1:3]
+# AAACAGCCAAGTGAAC.1 AAACAGCCACTGACCG.1 AAACAGCCATGATTGT.1
+#MIR1302.2HG 0 0 0
+#FAM138A 0 0 0
+#OR4F5 0 0 0
+#AL627309.1 0 0 0
+#AL627309.3 0 0 0
+#AL627309.2 0 0 0
+
+dim(MyRNAs)
+#[1] 36633 6326
+
+# make iCellR object
+my.obj <- make.obj(MyRNAs)
+
+# add ATAC-Seq data
+my.obj@atac.raw <- MyATAC
+my.obj@atac.main <- MyATAC
+
+# check your object
+my.obj
+
+
+###################################
+,--. ,-----. ,--.,--.,------.
+`--'' .--./ ,---. | || || .--. '
+,--.| | | .-. :| || || '--'.'
+| |' '--'\ --. | || || |
+`--' `-----' `----'`--'`--'`--' '--'
+###################################
+An object of class iCellR version: 1.6.2
+Raw/original data dimentions (rows,columns): 24127,6326
+Data conditions: no conditions/single sample
+Row names: MIR1302.2HG,TTLL10.AS1,MRPL20.AS1 ...
+Columns names: AAACAGCCAAGTGAAC.1,AAACAGCCACTGACCG.1,AAACAGCCATGATTGT.1 ...
+###################################
+ QC stats performed:FALSE, PCA performed:FALSE
+ Clustering performed:FALSE, Number of clusters:0
+ tSNE performed:FALSE, UMAP performed:FALSE, DiffMap performed:FALSE
+ Main data dimensions (rows,columns): 0,0
+ Normalization factors:,...
+ Imputed data dimensions (rows,columns):0,0
+############## scVDJ-seq ###########
+VDJ data dimentions (rows,columns):0,0
+############## CITE-seq ############
+ ADT raw data dimensions (rows,columns):0,0
+ ADT main data dimensions (rows,columns):0,0
+ ADT columns names:...
+ ADT row names:...
+############## scATAC-seq ############
+ ATAC raw data dimensions (rows,columns):21923,6326
+ ATAC main data dimensions (rows,columns):21923,6326
+ ATAC columns names:AAACAGCCAAGTGAAC.1...
+ ATAC row names:chr1.181218.181695...
+############## Spatial ###########
+Spatial data dimentions (rows,columns):0,0
+########### iCellR object ##########
+
+# from here do the regular scRNA-seq as expleind above
+
+# QC
+my.obj <- qc.stats(my.obj,
+ s.phase.genes = s.phase,
+ g2m.phase.genes = g2m.phase)
+
+# plot as mentioned above
+
+# filter
+my.obj <- cell.filter(my.obj,
+ min.mito = 0,
+ max.mito = 0.07 ,
+ min.genes = 500,
+ max.genes = 4000,
+ min.umis = 0,
+ max.umis = Inf)
+
+# normalize RNA
+my.obj <- norm.data(my.obj, norm.method = "ranked.glsf", top.rank = 500)
+
+# normalize ADT
+my.obj <- norm.adt(my.obj)
+
+# gene stats
+my.obj <- gene.stats(my.obj, which.data = "main.data")
+
+# find genes for PCA
+my.obj <- make.gene.model(my.obj, my.out.put = "data",
+ dispersion.limit = 1.5,
+ base.mean.rank = 500,
+ no.mito.model = T,
+ mark.mito = T,
+ interactive = F,
+ no.cell.cycle = T,
+ out.name = "gene.model")
+
+# run PCA and the rest is as above
+
+my.obj <- run.pca(my.obj, method = "gene.model", gene.list = my.obj@gene.model,data.type = "main")
+
+my.obj <- run.umap(my.obj, dims = 1:10)
+
+...
+ ```
+
+
+
+
+
+
+
+
+
+
+
From 125b1bb1557b1d5b584f7e73b6c47ef36419af20 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 11:00:36 -0400
Subject: [PATCH 09/74] Update README.md
---
README.md | 20 ++++++++++++++++++--
1 file changed, 18 insertions(+), 2 deletions(-)
diff --git a/README.md b/README.md
index 9dcd47e..3687911 100644
--- a/README.md
+++ b/README.md
@@ -2825,8 +2825,11 @@ VDJ data dimentions (rows,columns):0,0
Spatial data dimentions (rows,columns):0,0
########### iCellR object ##########
-# from here do the regular scRNA-seq as expleind above
+ ```
+From here do the regular scRNA-seq as expleind above. See example below
+
+ ```r
# QC
my.obj <- qc.stats(my.obj,
s.phase.genes = s.phase,
@@ -2866,9 +2869,22 @@ my.obj <- make.gene.model(my.obj, my.out.put = "data",
my.obj <- run.pca(my.obj, method = "gene.model", gene.list = my.obj@gene.model,data.type = "main")
+# tSNE
+my.obj <- run.pc.tsne(my.obj, dims = 1:10)
+
+# UMAP
my.obj <- run.umap(my.obj, dims = 1:10)
-...
+# KNetL
+my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap")
+
+# clustering based on KNetL
+
+my.obj <- iclust(my.obj, sensitivity = 150, data.type = "knetl")
+
+# clustering based on PCA
+
+# my.obj <- iclust(my.obj, sensitivity = 150, data.type = "pca", dims=1:10)
```
From 0b76efc036a8ddc9d230d2cd56ce41dad5a38140 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 11:44:06 -0400
Subject: [PATCH 10/74] Update README.md
---
README.md | 44 +++++++++++++++++++++++++++++++++++++++++---
1 file changed, 41 insertions(+), 3 deletions(-)
diff --git a/README.md b/README.md
index 3687911..65b7c4d 100644
--- a/README.md
+++ b/README.md
@@ -2876,15 +2876,53 @@ my.obj <- run.pc.tsne(my.obj, dims = 1:10)
my.obj <- run.umap(my.obj, dims = 1:10)
# KNetL
-my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap")
+my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 200, dim.redux = "umap")
# clustering based on KNetL
-my.obj <- iclust(my.obj, sensitivity = 150, data.type = "knetl")
+my.obj <- iclust(my.obj, sensitivity = 200, data.type = "knetl")
# clustering based on PCA
-# my.obj <- iclust(my.obj, sensitivity = 150, data.type = "pca", dims=1:10)
+# my.obj <- iclust(my.obj, sensitivity = 100, data.type = "pca", dims=1:10)
+
+# check clusters and adjust if needed (optinal)
+# cluster.plot(my.obj,plot.type = "knetl",interactive = F,cell.size = 0.5,cell.transparency = 1, anno.clust=T)
+# my.obj <- change.clust(my.obj, change.clust = 3, to.clust = 4)
+# my.obj <- change.clust(my.obj, change.clust = 3, to.clust = 10)
+
+# order clusters
+my.obj <- clust.ord(my.obj,top.rank = 500, how.to.order = "distance")
+
+
+# plot
+A= cluster.plot(my.obj,plot.type = "pca",interactive = F,cell.size = 0.5,cell.transparency = 1, anno.clust=T)
+B= cluster.plot(my.obj,plot.type = "umap",interactive = F,cell.size = 0.5,cell.transparency = 1, anno.clust=T)
+C= cluster.plot(my.obj,plot.type = "tsne",interactive = F,cell.size = 0.5,cell.transparency = 1, anno.clust=T)
+D= cluster.plot(my.obj,plot.type = "knetl",interactive = F,cell.size = 0.5,cell.transparency = 1, anno.clust=T)
+
+library(gridExtra)
+png('AllClusts.png', width = 12, height = 10, units = 'in', res = 300)
+grid.arrange(A,B,C,D)
+dev.off()
+
+# save object
+save(my.obj, file = "my.obj.Robj")
+
+# find markers
+marker.genes <- findMarkers(my.obj,
+ data.type = "main",
+ fold.change = 2,
+ padjval = 0.1,
+ uniq = F,
+ positive = T)
+
+MyGenes <- top.markers(marker.genes, topde = 10, min.base.mean = 0.2, filt.ambig = F)
+MyGenes <- unique(MyGenes)
+
+png('heatmap_gg_genes.png', width = 10, height = 10, units = 'in', res = 300)
+heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters",cell.sort = F, conds.to.plot = NULL)
+dev.off()
```
From 290947063d8643d5e242f4405dc9f0f9e7f45ea9 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 11:50:10 -0400
Subject: [PATCH 11/74] Update README.md
---
README.md | 7 +++++++
1 file changed, 7 insertions(+)
diff --git a/README.md b/README.md
index 65b7c4d..4243f73 100644
--- a/README.md
+++ b/README.md
@@ -2916,6 +2916,13 @@ marker.genes <- findMarkers(my.obj,
padjval = 0.1,
uniq = F,
positive = T)
+
+
+head(marker.genes)
+marker.genes1 <- cbind(row = rownames(marker.genes), marker.genes)
+head(marker.genes1)
+dim(marker.genes1)
+write.table((marker.genes1),file="marker.genes.tsv", sep="\t", row.names =F)
MyGenes <- top.markers(marker.genes, topde = 10, min.base.mean = 0.2, filt.ambig = F)
MyGenes <- unique(MyGenes)
From 5cdf6b93f1f48f1ab02c194e6f9b49ede5fd15c4 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 12:44:11 -0400
Subject: [PATCH 12/74] Update README.md
---
README.md | 21 +++++++++++++++++----
1 file changed, 17 insertions(+), 4 deletions(-)
diff --git a/README.md b/README.md
index 4243f73..ba80eca 100644
--- a/README.md
+++ b/README.md
@@ -2917,11 +2917,7 @@ marker.genes <- findMarkers(my.obj,
uniq = F,
positive = T)
-
-head(marker.genes)
marker.genes1 <- cbind(row = rownames(marker.genes), marker.genes)
-head(marker.genes1)
-dim(marker.genes1)
write.table((marker.genes1),file="marker.genes.tsv", sep="\t", row.names =F)
MyGenes <- top.markers(marker.genes, topde = 10, min.base.mean = 0.2, filt.ambig = F)
@@ -2932,7 +2928,24 @@ heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters"
dev.off()
```
+Work on scATAC data
+
+```r
+# normalize ACAT
+my.obj <- norm.data(my.obj, norm.method = "ranked.glsf", top.rank = 500, ATAC.data = TRUE, ATAC.filter = TRUE)
+marker.peaks <- findMarkers(my.obj,
+ data.type = "atac",
+ fold.change = 2,
+ padjval = 0.1,
+ uniq = F,
+ positive = T)
+
+marker.peaks1 <- cbind(row = rownames(marker.peaks), marker.genes)
+write.table((marker.peaks1),file="marker.peaks.tsv", sep="\t", row.names =F)
+
+head(marker.peaks1)
+```
From dad5cc944d725da6428dd9dd3f9184bffb3f8031 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 12:55:40 -0400
Subject: [PATCH 13/74] Update README.md
---
README.md | 52 +++++++++++++++++++++++++++++++++++++++++++++++++++-
1 file changed, 51 insertions(+), 1 deletion(-)
diff --git a/README.md b/README.md
index ba80eca..02466c5 100644
--- a/README.md
+++ b/README.md
@@ -2941,10 +2941,60 @@ marker.peaks <- findMarkers(my.obj,
uniq = F,
positive = T)
-marker.peaks1 <- cbind(row = rownames(marker.peaks), marker.genes)
+marker.peaks1 <- cbind(row = rownames(marker.peaks), marker.peaks)
write.table((marker.peaks1),file="marker.peaks.tsv", sep="\t", row.names =F)
head(marker.peaks1)
+# row baseMean baseSD
+#chr17.64986035.64986113 chr17.64986035.64986113 0.01217359 0.18257818
+#chr1.26542287.26542678 chr1.26542287.26542678 0.05828764 0.80077656
+#chr4.8199063.8199275 chr4.8199063.8199275 0.04280424 0.56205649
+#chr20.50274929.50275237 chr20.50274929.50275237 0.04684509 0.63361490
+#chr2.218382038.218382236 chr2.218382038.218382236 0.03122394 0.31153105
+#chr11.1760469.1760814 chr11.1760469.1760814 0.07050175 0.63322284
+# AvExpInCluster AvExpInOtherClusters foldChange
+#chr17.64986035.64986113 0.07868394 0.002346603 33.530999
+#chr1.26542287.26542678 0.33849093 0.016887273 20.044144
+#chr4.8199063.8199275 0.24803497 0.012481148 19.872769
+#chr20.50274929.50275237 0.27007513 0.013862584 19.482308
+#chr2.218382038.218382236 0.17736269 0.009631770 18.414340
+#chr11.1760469.1760814 0.39043782 0.023230813 16.806894
+# log2FoldChange pval padj clusters
+#chr17.64986035.64986113 5.067424 1.187697e-05 2.875415e-02 1
+#chr1.26542287.26542678 4.325109 3.653916e-05 8.539202e-02 1
+#chr4.8199063.8199275 4.312721 6.059691e-06 1.489472e-02 1
+#chr20.50274929.50275237 4.284093 2.871301e-05 6.779143e-02 1
+#chr2.218382038.218382236 4.202758 6.359572e-10 1.677019e-06 1
+#chr11.1760469.1760814 4.070981 6.813447e-11 1.800794e-07 1
+# gene cluster_1 cluster_2
+#chr17.64986035.64986113 chr17.64986035.64986113 0.0786839378 0.000000000
+#chr1.26542287.26542678 chr1.26542287.26542678 0.3384909326 0.008895062
+#chr4.8199063.8199275 chr4.8199063.8199275 0.2480349741 0.038672840
+#chr20.50274929.50275237 chr20.50274929.50275237 0.2700751295 0.028703704
+#chr2.218382038.218382236 chr2.218382038.218382236 0.1773626943 0.000000000
+#chr11.1760469.1760814 chr11.1760469.1760814 0.3904378238 0.004537037
+# cluster_3 cluster_4 cluster_5 cluster_6
+#chr17.64986035.64986113 0.000000000 0.0000000000 0.0006934750 0.0010038760
+#chr1.26542287.26542678 0.031485714 0.0029244992 0.0052261002 0.0041264535
+#chr4.8199063.8199275 0.000000000 0.0007226502 0.0027450683 0.0113212209
+#chr20.50274929.50275237 0.027092857 0.0121741140 0.0041820941 0.0051516473
+#chr2.218382038.218382236 0.004292857 0.0042095532 0.0006722307 0.0038561047
+#chr11.1760469.1760814 0.071678571 0.0177288136 0.0141820941 0.0061099806
+# cluster_7 cluster_8
+#chr17.64986035.64986113 0.003099029 0.009080357
+#chr1.26542287.26542678 0.047110680 0.044484375
+#chr4.8199063.8199275 0.014819417 0.028651786
+#chr20.50274929.50275237 0.025499029 0.028765625
+#chr2.218382038.218382236 0.011988350 0.035508929
+#chr11.1760469.1760814 0.042093204 0.050708705
+
+MyGenes <- top.markers(marker.peaks, topde = 10, min.base.mean = 0.2, filt.ambig = F)
+MyGenes <- unique(MyGenes)
+
+png('heatmap_gg_peaks.png', width = 10, height = 10, units = 'in', res = 300)
+heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters",cell.sort = F, conds.to.plot = NULL)
+dev.off()
+
```
From db3689aa9158d59a9aefcd766c6132dd730f8774 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 13:41:10 -0400
Subject: [PATCH 14/74] Update README.md
---
README.md | 8 +++++++-
1 file changed, 7 insertions(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 02466c5..4d3a8c7 100644
--- a/README.md
+++ b/README.md
@@ -2992,9 +2992,15 @@ MyGenes <- top.markers(marker.peaks, topde = 10, min.base.mean = 0.2, filt.ambig
MyGenes <- unique(MyGenes)
png('heatmap_gg_peaks.png', width = 10, height = 10, units = 'in', res = 300)
-heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters",cell.sort = F, conds.to.plot = NULL)
+heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters",cell.sort = F, conds.to.plot = NULL, data.type = "atac")
dev.off()
+my.obj <- run.impute(my.obj,data.type = "knetl", nn = 10, ATAC.data = FALSE)
+
+
+png('heatmap_gg_peaks.png', width = 10, height = 10, units = 'in', res = 300)
+heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters",cell.sort = F, conds.to.plot = NULL, data.type = "atac.imputed")
+dev.off()
```
From 9f46e1eb8a7bd4d03b34b508b6a76576ab5b3b38 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Fri, 16 Apr 2021 14:37:46 -0400
Subject: [PATCH 15/74] Update README.md
---
README.md | 61 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 61 insertions(+)
diff --git a/README.md b/README.md
index 4d3a8c7..db84961 100644
--- a/README.md
+++ b/README.md
@@ -3003,8 +3003,69 @@ heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters"
dev.off()
```
+Peak analysis
+```r
+# make bed file per cluster
+make.bed(marker.peaks)
+
+# load packages
+library(ChIPseeker)
+library(clusterProfiler)
+
+# load genome
+require(TxDb.Hsapiens.UCSC.hg38.knownGene)
+txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene
+Anno="org.Hs.eg.db"
+
+# load bed files
+Mylist1 = list.files(pattern=".bed")
+Mylist1
+
+Mylist <- as.list(Mylist1)
+NAMES <- gsub('_peaks.bed','',Mylist1)
+names(Mylist) <- NAMES
+files <- Mylist
+files
+
+# perform analysis (example)
+promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
+tagMatrixList <- lapply(files, getTagMatrix, windows=promoter)
+
+pdf("Plot_ProfileLineAll.pdf")
+plotAvgProf(tagMatrixList, xlim=c(-3000, 3000))
+dev.off()
+
+pdf('Plot_ProfileLine.pdf', width = 8, height = 10)
+plotAvgProf(tagMatrixList, xlim=c(-3000, 3000), facet="row")
+dev.off()
+
+pdf("Plot_heatmaps.pdf", width = 50, height = 6)
+tagHeatmap(tagMatrixList, xlim=c(-3000, 3000), color=NULL)
+dev.off()
+
+# annotate
+peakAnnoList <- lapply(files, annotatePeak, TxDb=txdb,
+ tssRegion=c(-3000, 3000), verbose=FALSE)
+# plot annotatin
+pdf("Plot_AnnoBar.pdf")
+plotAnnoBar(peakAnnoList)
+dev.off()
+
+############### peak annotation
+
+peakAnnoList <- lapply(files, annotatePeak, TxDb=txdb,
+ tssRegion=c(-3000, 3000), verbose=FALSE, annoDb=Anno)
+capture.output(peakAnnoList, file = "peakAnnoList.txt")
+
+genes = lapply(peakAnnoList, function(i) as.data.frame(i))
+
+lapply(1:length(genes), function(i) write.table(genes[[i]],
+ file = paste0(names(genes[i]), ".xls"),
+ row.names = FALSE, sep="\t"))
+
+```
From 398cba072da136674be98ec4ba5ebf12a10c97d8 Mon Sep 17 00:00:00 2001
From: Khodadadi-Jamayran
Date: Tue, 27 Apr 2021 12:47:32 -0400
Subject: [PATCH 16/74] scATAC
---
DESCRIPTION | 4 ++--
NAMESPACE | 1 +
R/F00100.R | 8 ++++++--
R/F0021.R | 8 +++++++-
R/F0022.R | 36 ++++++++++++++++++++++++++++-------
R/F0028.R | 11 +++++++++--
R/F0030.R | 8 +++++++-
R/F0031.R | 8 +++++++-
R/F0044.R | 1 +
R/F0057.R | 2 +-
R/F0058.R | 16 +++++++++++++++-
R/F0072.R | 26 +++++++++++++++++++++++++
R/F008.R | 43 ++++++++++++++++++++++++++++++++++++++----
man/clust.avg.exp.Rd | 2 +-
man/findMarkers.Rd | 2 +-
man/gene.plot.Rd | 2 +-
man/heatmap.gg.plot.Rd | 2 +-
man/hto.anno.Rd | 2 +-
man/make.bed.Rd | 17 +++++++++++++++++
man/norm.data.Rd | 14 +++++++++++++-
man/run.impute.Rd | 9 +++++++++
21 files changed, 194 insertions(+), 28 deletions(-)
create mode 100644 R/F0072.R
create mode 100644 man/make.bed.Rd
diff --git a/DESCRIPTION b/DESCRIPTION
index cbef2b6..d544417 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
Package: iCellR
Type: Package
Title: Analyzing High-Throughput Single Cell Sequencing Data
-Version: 1.6.2
+Version: 1.6.4
Authors@R: c(
person(given = 'Alireza', family = 'Khodadadi-Jamayran',role = c('aut','cre'), email = 'alireza.khodadadi.j@gmail.com', comment = c(ORCID = '0000-0003-2495-7504')),
person(given = 'Joseph', family = 'Pucella',role = c('aut','ctb'), comment = c(ORCID = '0000-0003-0875-8046')),
@@ -27,7 +27,7 @@ RoxygenNote: 7.1.1
BugReports: https://github.com/rezakj/iCellR/issues
URL: https://github.com/rezakj/iCellR
NeedsCompilation: yes
-Packaged: 2021-04-01 17:47:14 UTC; khodaa01
+Packaged: 2021-04-27 16:41:35 UTC; khodaa01
Author: Alireza Khodadadi-Jamayran [aut, cre]
(),
Joseph Pucella [aut, ctb] (),
diff --git a/NAMESPACE b/NAMESPACE
index 34a3847..3667903 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -36,6 +36,7 @@ export(iba)
export(iclust)
export(load.h5)
export(load10x)
+export(make.bed)
export(make.gene.model)
export(make.obj)
export(myImp)
diff --git a/R/F00100.R b/R/F00100.R
index 8136698..9c0ef8c 100644
--- a/R/F00100.R
+++ b/R/F00100.R
@@ -99,7 +99,11 @@ utils::globalVariables(c("%>%",
"wilcox.test",
"gsva",
"V5",
- "V6"))
+ "V6",
+ "AvExpInOtherClusters",
+ "X1",
+ "X2",
+ "X3"))
#########
# F001.load10x.R | F001.R
# F002.data.aggregation.R | F002.R
@@ -172,4 +176,4 @@ utils::globalVariables(c("%>%",
# add.10x.image.R | F0069.R
# spatial.plot.R | F0070.R
# i.score.R | F0071.R
-
+# make.bed.R | F0072.R
diff --git a/R/F0021.R b/R/F0021.R
index cbb1195..81c0f5d 100644
--- a/R/F0021.R
+++ b/R/F0021.R
@@ -2,7 +2,7 @@
#'
#' This function takes an object of class iCellR and creates an average gene expression for every cluster.
#' @param x An object of class iCellR.
-#' @param data.type Choose from "main" and "imputed", default = "main"
+#' @param data.type Choose from "main", "atac", "atac.imputed" and "imputed", default = "main"
#' @param conds.to.avg Choose the conditions you want to average, default = NULL (all conditions).
#' @param rounding.digits integer indicating the number of decimal places (round) or significant digits (signif) to be used.
#' @param round.num Rounding of Numbers, default = FALSE.
@@ -42,6 +42,12 @@ clust.avg.exp <- function (x = NULL,
if (data.type == "imputed") {
Table <- x@imputed.data
}
+ if (data.type == "atac") {
+ Table <- x@atac.main
+ }
+ if (data.type == "atac.imputed") {
+ Table <- x@atac.imputed
+ }
# Table = x@main.data
datalist <- list()
###
diff --git a/R/F0022.R b/R/F0022.R
index c2c2116..4d36ab0 100644
--- a/R/F0022.R
+++ b/R/F0022.R
@@ -3,6 +3,7 @@
#' This function takes an object of class iCellR and runs imputation on the main data.
#' @param x An object of class iCellR.
#' @param imp.method Choose between "iCellR.imp" and "magic", defualt = "iCellR.imp".
+#' @param ATAC.data If TURE, it would normalize ATAC-Seq data and not RNA-Seq, default = FALSE.
#' @param nn Number of neighboring cells to find, default = 10.
#' @param dims PC dimentions to be used for the analysis, default = 10.
#' @param data.type Choose between "tsne", "pca", "umap", "diffusion", "knetl", default = "pca".
@@ -17,22 +18,35 @@
#' @param verbose 'int' or 'boolean', optional (default : 1) If 'TRUE' or '> 0', message verbose updates.
#' @param n.jobs 'int', optional (default: 1) The number of jobs to use for the computation. If -1 all CPUs are used. If 1 is given, no parallel computing code is used at all, which is useful for debugging. For n_jobs below -1, (n.cpus + 1 + n.jobs) are used. Thus for n_jobs = -2, all CPUs but one are used
#' @param seed int or 'NULL', random state (default: 'NULL')
+#' @param rounding.digits integer indicating the number of decimal places (round) or significant digits (signif) to be used.
+#' @param round.num Rounding of Numbers, default = FALSE.
#' @return An object of class iCellR.
#' @import progress
#' @export
run.impute <- function (x = NULL,
- imp.method = "iCellR.imp", dims = 1:10, nn = 10,
- data.type = "pca",genes = "all_genes", k = 10, alpha = 15, t = "auto",
+ imp.method = "iCellR.imp",
+ dims = 1:10, nn = 10, ATAC.data = FALSE,
+ rounding.digits = 4,
+ round.num = TRUE,
+ data.type = "pca",genes = "all_genes",
+ k = 10, alpha = 15, t = "auto",
npca = 100, init = NULL, t.max = 20,
- knn.dist.method = "euclidean", verbose = 1, n.jobs = 1,
+ knn.dist.method = "euclidean",
+ verbose = 1, n.jobs = 1,
seed = NULL) {
if ("iCellR" != class(x)[1]) {
stop("x should be an object of class iCellR")
}
# get data
start_time1 <- Sys.time()
- #####
+ ####
+ # get data RNA
DATA <- x@main.data
+ # get data ATAC
+ if (ATAC.data == TRUE) {
+ DATA <- x@atac.main
+ }
+ ####
message(paste(" main data dimentions:",dim(DATA)[1],"genes and",dim(DATA)[2],"cells"))
message(" Genes with no coverage are being removed from the matrix ...")
DATA <- DATA[ rowSums(DATA) > 0, ]
@@ -142,9 +156,17 @@ run.impute <- function (x = NULL,
knn.dist.method = knn.dist.method, verbose = verbose, n.jobs = n.jobs,
seed = seed)
DATA <- as.data.frame(t(data_MAGIC$result))
- # return
- DATA <- round(DATA, digits = 3)
- attributes(x)$imputed.data <- DATA
+ #####
+ if (round.num == TRUE) {
+ DATA <- round(DATA, digits = rounding.digits)
+ }
+#########
+ if (ATAC.data == TRUE) {
+ attributes(x)$atac.imputed <- DATA
+ }
+ if (ATAC.data == FALSE) {
+ attributes(x)$imputed.data <- DATA
+ }
# return
return(x)
}
diff --git a/R/F0028.R b/R/F0028.R
index 08237f6..0278ed1 100644
--- a/R/F0028.R
+++ b/R/F0028.R
@@ -2,7 +2,7 @@
#'
#' This function takes an object of class iCellR and performs differential expression (DE) analysis to find marker genes for each cluster.
#' @param x An object of class iCellR.
-#' @param data.type Choose from "main" and "imputed", default = "main"
+#' @param data.type Choose from "main", "atac", "atac.imputed" and "imputed", default = "main"
#' @param fold.change A number that designates the minimum fold change for out put, default = 2.
#' @param pval.test Choose from "t.test", "wilcox.test", default = "t.test".
#' @param p.adjust.method Correction method. Choose from "holm", "hochberg", "hommel", "bonferroni", "BH", "BY","fdr", "none", default = "hochberg".
@@ -34,9 +34,16 @@ findMarkers <- function (x = NULL,
if (data.type == "imputed") {
dat <- x@imputed.data
}
+ if (data.type == "atac") {
+ dat <- x@atac.main
+ }
+ if (data.type == "atac.imputed") {
+ dat <- x@atac.imputed
+ }
# get cluster data
# get avrages
- x <- clust.avg.exp(x, data.type = data.type)
+ x <- clust.avg.exp(x, data.type = data.type)
+##########
DATA <- x@best.clust
if(!is.numeric(DATA$clusters)){
stop("Cluster names have to be numeric")
diff --git a/R/F0030.R b/R/F0030.R
index 8978200..d14c398 100644
--- a/R/F0030.R
+++ b/R/F0030.R
@@ -10,7 +10,7 @@
#' @param interactive If TRUE an html interactive file will be made, default = TRUE.
#' @param out.name Output name for html file if interactive = TRUE, default = "plot".
#' @param no.key If you want a color legend key, default = FALSE.
-#' @param data.type Choose from "main" and "imputed", default = "main".
+#' @param data.type Choose from "main", "atac", atac.imputed and "imputed", default = "main".
#' @param min.scale Set a minimum color scale, default = -2.5.
#' @param max.scale Set a maximum color scale, default = 2.5.
#' @param cex.col Chhose a size, default = 10.
@@ -61,6 +61,12 @@ heatmap.gg.plot <- function (x = NULL,
if (data.type == "imputed") {
DATAmain <- x@imputed.data
}
+ if (data.type == "atac") {
+ DATAmain <- x@atac.main
+ }
+ if (data.type == "atac.imputed") {
+ DATAmain <- x@atac.imputed
+ }
AllGenes = row.names(DATAmain)
absent = which((gene %in% AllGenes) == FALSE)
absentgenes = gene[absent]
diff --git a/R/F0031.R b/R/F0031.R
index a594fa9..68cecd2 100644
--- a/R/F0031.R
+++ b/R/F0031.R
@@ -8,7 +8,7 @@
#' @param plot.data.type Choose between "tsne", "pca", "umap", "knetl", "diffusion", "pseudo.A" and "pseudo.B", default = "tsne".
#' @param clust.dim 2 for 2D plots and 3 for 3D plots, default = 2.
#' @param col.by Choose from "clusters" and "conditions", default = "clusters".
-#' @param data.type Choose from "main" or "imputed", default = "main".
+#' @param data.type Choose from "main", "atac, "atac.imputed" and "imputed", default = "main".
#' @param scaleValue Scale the colors, default = FALSE.
#' @param min.scale If scaleValue = TRUE, set a number for min, default = -2.5.
#' @param max.scale If scaleValue = TRUE, set a number for max, default = 2.5.
@@ -88,6 +88,12 @@ gene.plot <- function (x = NULL,
if (data.type == "imputed") {
DATAmain <- x@imputed.data
}
+ if (data.type == "atac") {
+ DATAmain <- x@atac.main
+ }
+ if (data.type == "atac.imputed") {
+ DATAmain <- x@atac.imputed
+ }
##########
#######
if (is.null(gene)) {
diff --git a/R/F0044.R b/R/F0044.R
index 3f27fc4..2c26641 100644
--- a/R/F0044.R
+++ b/R/F0044.R
@@ -6,6 +6,7 @@ setClass("iCellR", representation (raw.data = "data.frame",
metadata = "data.frame",
atac.raw = "data.frame",
atac.main = "data.frame",
+ atac.imputed = "data.frame",
spatial.info = "list",
spatial.data = "data.frame",
stats = "data.frame",
diff --git a/R/F0057.R b/R/F0057.R
index cf7e6b5..211848d 100644
--- a/R/F0057.R
+++ b/R/F0057.R
@@ -2,7 +2,7 @@
#'
#' Demultiplexing HTOs
#' @param hto.data HTO raw data
-#' @param cov.thr A number which avrage coverage is devided by to set a thershold for low coverage, default = 10.
+#' @param cov.thr A number which average coverage is divided by to set a threshold for low coverage. For example 10 means it is 10 time less than the average. default = 10.
#' @param assignment.thr A percent above which you decide to set as a good sample assignment/HTO, default = 80.
#' @return An object of class iCellR
#' @examples
diff --git a/R/F0058.R b/R/F0058.R
index f90412b..19e8a68 100644
--- a/R/F0058.R
+++ b/R/F0058.R
@@ -35,7 +35,21 @@ run.knetl <- function (x = NULL,
do.redux = TRUE,
run.iclust = FALSE,
return.graph = FALSE) {
- #
+ #####
+ if ("iCellR" != class(x)[1]) {
+ stop("x should be an object of class iCellR")
+ }
+ #####
+ message("################# IMPORTANT NOTE ##################","")
+ message(paste(" Zoom is set to:", zoom,"! Make sure it is good for this data."))
+ message("
+# For data with less than 5000 cells use a zoom of about 100-200.
+# For data with 5000-10000 cells use a zoom of about 100-300.
+# For data with 10000-30000 cells use a zoom of about 200-400.
+# For data with more than 30000 cells use a zoom of about 400-600.
+# A zoom value of 400 and dims = 1:20 is usually good for big data but adjust it for intended resolution.")
+ message("###################################################","")
+ ######
start_time1 <- Sys.time()
# cluster
if(data.type == "pca") {
diff --git a/R/F0072.R b/R/F0072.R
new file mode 100644
index 0000000..ebbc141
--- /dev/null
+++ b/R/F0072.R
@@ -0,0 +1,26 @@
+#' Make BED Files
+#'
+#' This function takes peak marker files and makes the bed files per cluster.
+#' @param x Peak marker file.
+#' @return Bed files
+#' @export
+make.bed <- function (x = NULL) {
+# get filed
+ results <- subset(x,select=c(gene,clusters,AvExpInOtherClusters))
+ peaks <- as.character(results$gene)
+ peaks <- (gsub("\\.","_",peaks))
+ peaks <- data.frame(do.call('rbind', strsplit(as.character(peaks),'_',fixed=TRUE)))
+ results <- cbind(peaks,results)
+#######
+ My.clusters <- unique(results$clusters)
+#######
+ for(i in My.clusters){
+ dat <- subset(results, results$clusters == i)
+ dat <- subset(dat,select=c(X1,X2,X3,gene,AvExpInOtherClusters))
+ MyName <- paste("peaks_cluster_",i,".bed",sep="")
+ if(dim(dat)[1] > 0) {
+ write.table(dat, MyName, sep="\t", quote = FALSE,row.names =FALSE, col.names = FALSE,)
+ }
+ }
+}
+
diff --git a/R/F008.R b/R/F008.R
index f5956fe..c5e2125 100644
--- a/R/F008.R
+++ b/R/F008.R
@@ -7,6 +7,10 @@
#' @param top.rank If the method is set to "ranked.glsf", you need to set top number of genes sorted based on global base mean, default = 500.
#' @param spike.in.factors A numeric vector of spike-in values with the same cell id order as the main data.
#' @param rpm.factor If the norm.method is set to "rpm" the library sizes would be divided by this number, default = 1000 (higher numbers recomanded for bulk RNA-Seq).
+#' @param rounding.digits integer indicating the number of decimal places (round) or significant digits (signif) to be used.
+#' @param round.num Rounding of Numbers, default = FALSE.
+#' @param ATAC.data If TURE, it would normalize ATAC-Seq data and not RNA-Seq, default = FALSE.
+#' @param ATAC.filter If TURE, all the cells filtered in RNA-Seq will be filtered in ATAC-Seq. This needs to be done for both data to match, default = TRUE.
#' @return An object of class iCellR.
#' @examples
#'
@@ -17,11 +21,22 @@ norm.data <- function (x = NULL,
norm.method = "ranked.glsf",
top.rank = 500,
spike.in.factors = NULL,
- rpm.factor = 1000) {
+ rpm.factor = 1000,
+ rounding.digits = 3,
+ round.num = TRUE,
+ ATAC.data = FALSE,
+ ATAC.filter = TRUE) {
if ("iCellR" != class(x)[1]) {
stop("x should be an object of class iCellR")
}
+#######
+# get data RNA
DATA <- x@main.data
+# get data ATAC
+ if (ATAC.data == TRUE) {
+ DATA <- x@atac.main
+ }
+#######
if (norm.method == "global.glsf") {
libSiz <- colSums(DATA)
norm.facts <- as.numeric(libSiz) / mean(as.numeric(libSiz))
@@ -56,9 +71,29 @@ norm.data <- function (x = NULL,
normalized <- as.data.frame(sweep(dataMat, 2, norm.facts, `/`))
}
#####
- normalized <- round(normalized, digits = 3)
- attributes(x)$main.data <- normalized
- attributes(x)$norm.factors <- norm.facts
+ if (round.num == TRUE) {
+ normalized <- round(normalized, digits = rounding.digits)
+ }
+######
+ # get data ATAC
+ if (ATAC.data == TRUE) {
+ if (ATAC.filter == TRUE) {
+ dat = normalized
+ MyIDsToKeep <- colnames(x@main.data)
+ normalized <- dat[ , which(names(dat) %in% MyIDsToKeep)]
+ #### order colnames by the original data
+ normalized <- as.matrix(normalized)
+ normalized <- normalized[ , order(match(colnames(normalized),MyIDsToKeep))]
+ normalized <- as.data.frame(normalized)
+ }
+ attributes(x)$atac.main <- normalized
+ }
+ # get data ATAC
+ if (ATAC.data == FALSE) {
+ attributes(x)$main.data <- normalized
+ attributes(x)$norm.factors <- norm.facts
+ }
+############
return(x)
}
diff --git a/man/clust.avg.exp.Rd b/man/clust.avg.exp.Rd
index 97c58c5..3635aa1 100644
--- a/man/clust.avg.exp.Rd
+++ b/man/clust.avg.exp.Rd
@@ -15,7 +15,7 @@ clust.avg.exp(
\arguments{
\item{x}{An object of class iCellR.}
-\item{data.type}{Choose from "main" and "imputed", default = "main"}
+\item{data.type}{Choose from "main", "atac", "atac.imputed" and "imputed", default = "main"}
\item{conds.to.avg}{Choose the conditions you want to average, default = NULL (all conditions).}
diff --git a/man/findMarkers.Rd b/man/findMarkers.Rd
index 03df140..df46fd3 100644
--- a/man/findMarkers.Rd
+++ b/man/findMarkers.Rd
@@ -19,7 +19,7 @@ findMarkers(
\arguments{
\item{x}{An object of class iCellR.}
-\item{data.type}{Choose from "main" and "imputed", default = "main"}
+\item{data.type}{Choose from "main", "atac", "atac.imputed" and "imputed", default = "main"}
\item{pval.test}{Choose from "t.test", "wilcox.test", default = "t.test".}
diff --git a/man/gene.plot.Rd b/man/gene.plot.Rd
index 62ae690..eba5305 100644
--- a/man/gene.plot.Rd
+++ b/man/gene.plot.Rd
@@ -41,7 +41,7 @@ gene.plot(
\item{conds.to.plot}{Choose the conditions you want to see in the plot, default = NULL (all conditions).}
-\item{data.type}{Choose from "main" or "imputed", default = "main".}
+\item{data.type}{Choose from "main", "atac, "atac.imputed" and "imputed", default = "main".}
\item{box.to.test}{A cluster number so that all the boxes in the box plot would be compared to. If set to "0" the cluster with the highest avrage would be choosen, default = 0.}
diff --git a/man/heatmap.gg.plot.Rd b/man/heatmap.gg.plot.Rd
index 8e7cc29..6113b69 100644
--- a/man/heatmap.gg.plot.Rd
+++ b/man/heatmap.gg.plot.Rd
@@ -28,7 +28,7 @@ heatmap.gg.plot(
\item{cell.sort}{If FALSE the cells will not be sorted based on their distance, default = TRUE.}
-\item{data.type}{Choose from "main" and "imputed", default = "main".}
+\item{data.type}{Choose from "main", "atac", atac.imputed and "imputed", default = "main".}
\item{cluster.by}{Choose from "clusters" or "none", default = "clusters".}
diff --git a/man/hto.anno.Rd b/man/hto.anno.Rd
index 703221a..67ffc8b 100644
--- a/man/hto.anno.Rd
+++ b/man/hto.anno.Rd
@@ -9,7 +9,7 @@ hto.anno(hto.data = "data.frame", cov.thr = 10, assignment.thr = 80)
\arguments{
\item{hto.data}{HTO raw data}
-\item{cov.thr}{A number which avrage coverage is devided by to set a thershold for low coverage, default = 10.}
+\item{cov.thr}{A number which average coverage is divided by to set a threshold for low coverage. For example 10 means it is 10 time less than the average. default = 10.}
\item{assignment.thr}{A percent above which you decide to set as a good sample assignment/HTO, default = 80.}
}
diff --git a/man/make.bed.Rd b/man/make.bed.Rd
new file mode 100644
index 0000000..e927886
--- /dev/null
+++ b/man/make.bed.Rd
@@ -0,0 +1,17 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/F0072.R
+\name{make.bed}
+\alias{make.bed}
+\title{Make BED Files}
+\usage{
+make.bed(x = NULL)
+}
+\arguments{
+\item{x}{Peak marker file.}
+}
+\value{
+Bed files
+}
+\description{
+This function takes peak marker files and makes the bed files per cluster.
+}
diff --git a/man/norm.data.Rd b/man/norm.data.Rd
index 70a12af..aa9d31e 100644
--- a/man/norm.data.Rd
+++ b/man/norm.data.Rd
@@ -9,7 +9,11 @@ norm.data(
norm.method = "ranked.glsf",
top.rank = 500,
spike.in.factors = NULL,
- rpm.factor = 1000
+ rpm.factor = 1000,
+ rounding.digits = 3,
+ round.num = TRUE,
+ ATAC.data = FALSE,
+ ATAC.filter = TRUE
)
}
\arguments{
@@ -23,6 +27,14 @@ Choose from "global.glsf", "ranked.glsf","spike.in" or no.norm, default = "ranke
\item{spike.in.factors}{A numeric vector of spike-in values with the same cell id order as the main data.}
\item{rpm.factor}{If the norm.method is set to "rpm" the library sizes would be divided by this number, default = 1000 (higher numbers recomanded for bulk RNA-Seq).}
+
+\item{rounding.digits}{integer indicating the number of decimal places (round) or significant digits (signif) to be used.}
+
+\item{round.num}{Rounding of Numbers, default = FALSE.}
+
+\item{ATAC.data}{If TURE, it would normalize ATAC-Seq data and not RNA-Seq, default = FALSE.}
+
+\item{ATAC.filter}{If TURE, all the cells filtered in RNA-Seq will be filtered in ATAC-Seq. This needs to be done for both data to match, default = TRUE.}
}
\value{
An object of class iCellR.
diff --git a/man/run.impute.Rd b/man/run.impute.Rd
index 94a4b45..5d38c51 100644
--- a/man/run.impute.Rd
+++ b/man/run.impute.Rd
@@ -9,6 +9,9 @@ run.impute(
imp.method = "iCellR.imp",
dims = 1:10,
nn = 10,
+ ATAC.data = FALSE,
+ rounding.digits = 4,
+ round.num = TRUE,
data.type = "pca",
genes = "all_genes",
k = 10,
@@ -32,6 +35,12 @@ run.impute(
\item{nn}{Number of neighboring cells to find, default = 10.}
+\item{ATAC.data}{If TURE, it would normalize ATAC-Seq data and not RNA-Seq, default = FALSE.}
+
+\item{rounding.digits}{integer indicating the number of decimal places (round) or significant digits (signif) to be used.}
+
+\item{round.num}{Rounding of Numbers, default = FALSE.}
+
\item{data.type}{Choose between "tsne", "pca", "umap", "diffusion", "knetl", default = "pca".}
\item{genes}{character or integer vector, default: NULL vector of column names or column indices for which to return smoothed data If 'all_genes' or NULL, the entire smoothed matrix is returned}
From 4cf6119dab0dab08c6000039ee097a43c5f89a04 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 27 Apr 2021 14:28:57 -0400
Subject: [PATCH 17/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index db84961..7cf021b 100644
--- a/README.md
+++ b/README.md
@@ -4,7 +4,7 @@
[![Build Status](https://travis-ci.com/rezakj/iCellR.svg?branch=master)](https://travis-ci.com/rezakj/iCellR)
# iCellR
-iCellR is an interactive R package to work with high-throughput single cell sequencing technologies (i.e scRNA-seq, scVDJ-seq, ST and CITE-seq).
+iCellR is an interactive R package to work with high-throughput single cell sequencing technologies (i.e scRNA-seq, scVDJ-seq, scATAC-seq, CITE-Seq and Spatial Transcriptomics (ST)).
### News (July 2020): See iCellR version 1.5.5 with new cell cycle analysis for G0, G1S, G2M, M, G1M and S [phase](https://genome.med.nyu.edu/results/external/iCellR/example1/All_cellcycle.png), Pseudotime Abstract KNetL map [(PAK map)](https://genome.med.nyu.edu/results/external/iCellR/example1/pseudotime.KNetL.png) and gene-gene [correlations](https://genome.med.nyu.edu/results/external/iCellR/example1/gene-gene.correlation.png). See below for how to.
From 24e81f199da3ad63f12f2f5b91fde3f5f46e11a0 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 27 Apr 2021 14:33:49 -0400
Subject: [PATCH 18/74] Update README.md
---
README.md | 2 ++
1 file changed, 2 insertions(+)
diff --git a/README.md b/README.md
index 7cf021b..1ce28b2 100644
--- a/README.md
+++ b/README.md
@@ -6,6 +6,8 @@
# iCellR
iCellR is an interactive R package to work with high-throughput single cell sequencing technologies (i.e scRNA-seq, scVDJ-seq, scATAC-seq, CITE-Seq and Spatial Transcriptomics (ST)).
+### News (April 2021): Use iCellR version 1.6.4 for scATAC-seq and Spatial Transcriptomics (ST). Use function i.score for scoring (scoring cells based on genes of your interest) methods (i.e. tirosh, mean, sum, gsva, ssgsea, zscore and plage).
+
### News (July 2020): See iCellR version 1.5.5 with new cell cycle analysis for G0, G1S, G2M, M, G1M and S [phase](https://genome.med.nyu.edu/results/external/iCellR/example1/All_cellcycle.png), Pseudotime Abstract KNetL map [(PAK map)](https://genome.med.nyu.edu/results/external/iCellR/example1/pseudotime.KNetL.png) and gene-gene [correlations](https://genome.med.nyu.edu/results/external/iCellR/example1/gene-gene.correlation.png). See below for how to.
### News (May 2020): see our dimensionality reduction called [KNetL map](https://genome.med.nyu.edu/results/external/iCellR/example1/Allclusts.Annotated.png) (pronounced like "nettle"). [KNetL](https://www.biorxiv.org/content/10.1101/2020.05.05.078550v1.full) map is capable of zooming and shows a lot more details compared to tSNE and UMAP.
From 4c39cafdb984ea24e78995d5c81963cef899d714 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 27 Apr 2021 15:03:37 -0400
Subject: [PATCH 19/74] Update README.md
---
README.md | 13 +++++++++----
1 file changed, 9 insertions(+), 4 deletions(-)
diff --git a/README.md b/README.md
index 1ce28b2..2c6d620 100644
--- a/README.md
+++ b/README.md
@@ -429,13 +429,18 @@ my.obj <- run.umap(my.obj, dims = 1:10)
my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap") # (Important note!) don't forget to set the zoom in the right range
########################### IMPORTANT NOTE ########################################
-#### Because KNetl has a very high resolution it's best to use a dim of 20 (this usually works best for most data)
#### For zooming use the zoom option.
-# For data with less than 5000 cells use a zoom of about 100-200.
+# For data with less than 5000 cells use a zoom of about 100-200.
# For data with 5000-10000 cells use a zoom of about 100-300.
-# For data with 10000-30000 cells use a zoom of about 200-400.
+# For data with 10000-30000 cells use a zoom of about 200-500.
# For data with more than 30000 cells use a zoom of about 400-600.
-#### A zoom value of 400 is usually good for big data but adjust it for intended resolution.
+# zoom 400 is usually good for big data but adjust for intended resolution.
+# Lower number for zoom in higher for zoom out (its reverse).
+# dims = 1:20 is generally good for most data.
+# other parameters are best as default.
+
+ ***KNetL map is very dynamic with zoom and dims***
+
#### Just like a microscope, you need to zoom to see the intended amount of details.
#### Here we use a zoom of 100 or 110 but this might not be ideal for your data.
#### example: # my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 400)
From c48e914c5645cbf953f49b8da867ff61189eae2b Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Thu, 20 May 2021 12:07:39 -0400
Subject: [PATCH 20/74] Update README.md
---
README.md | 12 ++++++------
1 file changed, 6 insertions(+), 6 deletions(-)
diff --git a/README.md b/README.md
index 2c6d620..f4ee684 100644
--- a/README.md
+++ b/README.md
@@ -428,19 +428,19 @@ my.obj <- run.umap(my.obj, dims = 1:10)
# Because knetl has a very high resolution it's best to use a dim of 20 (this usually works best for most data)
my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 110, dim.redux = "umap") # (Important note!) don't forget to set the zoom in the right range
-########################### IMPORTANT NOTE ########################################
-#### For zooming use the zoom option.
-# For data with less than 5000 cells use a zoom of about 100-200.
+########################## IMPORTANT DISCLAIMER NOTE ###########################
+ *** KNetL map is very dynamic with zoom and dims! ***
+ *** Therefore it needs to be adjusted! ***
+# For data with less than 1000 cells use a zoom of about 5-50.
+# For data with 1000-5000 cells use a zoom of about 50-200.
# For data with 5000-10000 cells use a zoom of about 100-300.
# For data with 10000-30000 cells use a zoom of about 200-500.
# For data with more than 30000 cells use a zoom of about 400-600.
# zoom 400 is usually good for big data but adjust for intended resolution.
-# Lower number for zoom in higher for zoom out (its reverse).
+# Lower number for zoom in and higher for zoom out (its reverse).
# dims = 1:20 is generally good for most data.
# other parameters are best as default.
- ***KNetL map is very dynamic with zoom and dims***
-
#### Just like a microscope, you need to zoom to see the intended amount of details.
#### Here we use a zoom of 100 or 110 but this might not be ideal for your data.
#### example: # my.obj <- run.knetl(my.obj, dims = 1:20, zoom = 400)
From 26706dc85a415401ffe51f076eaf8ccc4d463ba4 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:01:18 -0400
Subject: [PATCH 21/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index f4ee684..85aa2a1 100644
--- a/README.md
+++ b/README.md
@@ -508,7 +508,7 @@ In this function we provide a variety of many other options for you to explore t
| ------------- | ------------- | ------------- |
| ward.D, ward.D2, single, complete, average, mcquitty, median, centroid, kmeans| euclidean, maximum, manhattan, canberra, binary, minkowski or NULL | kl, ch, hartigan, ccc, scott, marriot, trcovw, tracew, friedman, rubin, cindex, db, silhouette, duda, pseudot2, beale, ratkowsky, ball, ptbiserial, gap, frey, mcclain, gamma, gplus, tau, dunn, hubert, sdindex, dindex, sdbw |
-### Conventionally people cluster based on PCA data however because KNetL map is more powerful we recommend clustering based on KNetL map.
+### Conventionally people cluster based on PCA data (usually first 10 dimensions) however you have the option of choosing tSNA, UMAP and KNetL map as well. If you have adjusted your KNetL map and are confident about the results we recommend clustering based on KNetL map.
This is one of the harder parts of the analysis and sometimes you need to adjust your clustering based on marker genes. This means you might need to merge some clusters, gate (see our cell gating tools) or try different sensitivities to find more or less communities.
From 84d5ddcc2defa7c8c415cda148bff0b0ea21853f Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:14:30 -0400
Subject: [PATCH 22/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 85aa2a1..f0c0fe2 100644
--- a/README.md
+++ b/README.md
@@ -113,7 +113,7 @@ To see the help page for each function use question mark as:
- Aggregate data
-Conditions in iCellR are set in the header of the data and are separated by an underscore (_). Let's say you want to merge multiple datasets and run iCellR in aggregate mode. Here’s an example: I divided this sample into three sets and then aggregate them into one matrix.
+## Conditions in iCellR are set or shown in the column names of the data and are separated by an underscore (_). Let's say you want to merge multiple datasets (data frames/matrices) into one file and run iCellR in aggregate mode (all samples together). You can do so using "data.aggregation" function. Here’s an example: I divided this sample into four datasets and then aggregated them into one matrix. Here we are assuming you have four samples (e.g. WT,KO,Ctrl,KD). In this way, iCellR will know you have 4 samples for the rest of the analysis (e.g. batch alignment, plots, DE, etc.).
```r
dim(my.data)
From fdb8d103fe62c391b19cdd3fc059ad48751c0e8c Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:15:07 -0400
Subject: [PATCH 23/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index f0c0fe2..afed3b2 100644
--- a/README.md
+++ b/README.md
@@ -113,7 +113,7 @@ To see the help page for each function use question mark as:
- Aggregate data
-## Conditions in iCellR are set or shown in the column names of the data and are separated by an underscore (_). Let's say you want to merge multiple datasets (data frames/matrices) into one file and run iCellR in aggregate mode (all samples together). You can do so using "data.aggregation" function. Here’s an example: I divided this sample into four datasets and then aggregated them into one matrix. Here we are assuming you have four samples (e.g. WT,KO,Ctrl,KD). In this way, iCellR will know you have 4 samples for the rest of the analysis (e.g. batch alignment, plots, DE, etc.).
+## Conditions in iCellR are set or shown in the column names of the data and are separated by an underscore ( "_" ). Let's say you want to merge multiple datasets (data frames/matrices) into one file and run iCellR in aggregate mode (all samples together). You can do so using "data.aggregation" function. Here’s an example: I divided this sample into four datasets and then aggregated them into one matrix. Here we are assuming you have four samples (e.g. WT,KO,Ctrl,KD). In this way, iCellR will know you have 4 samples for the rest of the analysis (e.g. batch alignment, plots, DE, etc.).
```r
dim(my.data)
From 69c094e941918360431c7d4a68e199118fe19bcf Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:15:44 -0400
Subject: [PATCH 24/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index afed3b2..84252c1 100644
--- a/README.md
+++ b/README.md
@@ -113,7 +113,7 @@ To see the help page for each function use question mark as:
- Aggregate data
-## Conditions in iCellR are set or shown in the column names of the data and are separated by an underscore ( "_" ). Let's say you want to merge multiple datasets (data frames/matrices) into one file and run iCellR in aggregate mode (all samples together). You can do so using "data.aggregation" function. Here’s an example: I divided this sample into four datasets and then aggregated them into one matrix. Here we are assuming you have four samples (e.g. WT,KO,Ctrl,KD). In this way, iCellR will know you have 4 samples for the rest of the analysis (e.g. batch alignment, plots, DE, etc.).
+## Conditions in iCellR are set or shown in the column names of the data and are separated by an underscore "_" sign. Let's say you want to merge multiple datasets (data frames/matrices) into one file and run iCellR in aggregate mode (all samples together). You can do so using "data.aggregation" function. Here’s an example: I divided this sample into four datasets and then aggregated them into one matrix. Here we are assuming you have four samples (e.g. WT,KO,Ctrl,KD). In this way, iCellR will know you have 4 samples for the rest of the analysis (e.g. batch alignment, plots, DE, etc.).
```r
dim(my.data)
From bfeef17ad18e824ebf7148ec487915aa3cf9b0ca Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:19:40 -0400
Subject: [PATCH 25/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 84252c1..9af7195 100644
--- a/README.md
+++ b/README.md
@@ -271,7 +271,7 @@ This step is optional and is for having the same number of cells for each condit
- Normalize data
-You have a few options to normalize your data based on your study. You can also normalize your data using tools other than iCellR and import your data to iCellR. We recommend "ranked.glsf" normalization for most single cell studies. This normalization is great for fixing matrixes with lots of zeros and because it's geometric it is great for fixing for batch effects, as long as all the data is aggregated into one file (to aggregate your data see "aggregating data" section above).
+You have a few options to normalize your data based on your study. You can also normalize your data using tools other than iCellR and import your data to iCellR. We recommend "ranked.glsf" normalization for most single cell studies. This normalization is great for fixing matrixes with lots of zeros and because it's geometric it will reduce some of batch differences in the library sizes, as long as all the data is aggregated into one file (to aggregate your data see "aggregating data" section above).
```r
my.obj <- norm.data(my.obj,
From d79089b202e107792d423ce0d386fefa4c658366 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:24:40 -0400
Subject: [PATCH 26/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 9af7195..eb66d3d 100644
--- a/README.md
+++ b/README.md
@@ -271,7 +271,7 @@ This step is optional and is for having the same number of cells for each condit
- Normalize data
-You have a few options to normalize your data based on your study. You can also normalize your data using tools other than iCellR and import your data to iCellR. We recommend "ranked.glsf" normalization for most single cell studies. This normalization is great for fixing matrixes with lots of zeros and because it's geometric it will reduce some of batch differences in the library sizes, as long as all the data is aggregated into one file (to aggregate your data see "aggregating data" section above).
+## You have a few options to normalize your data based on your study. You can also normalize your data using tools other than iCellR and import your data to iCellR. We recommend "ranked.glsf" normalization for most single cell studies. This normalization is great for fixing matrixes with lots of zeros and because it's geometric it will reduce some of batch differences in the library sizes, as long as all the data is aggregated into one file (to aggregate your data see "aggregating data" section above). GLSF stands for Geometric Library Size Factor, this is very similar to the normalization done by (DESeq2)[https://bioconductor.org/packages/release/bioc/html/DESeq2.html] and the ranked part would take the sum of the top most expressed genes as your library size instead of the full LB size which is to help resduce some of the drop out effects on normalization.
```r
my.obj <- norm.data(my.obj,
From eb149a8891a8fba2c044a440d24775cae1d604f4 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:25:14 -0400
Subject: [PATCH 27/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index eb66d3d..f027fe1 100644
--- a/README.md
+++ b/README.md
@@ -271,7 +271,7 @@ This step is optional and is for having the same number of cells for each condit
- Normalize data
-## You have a few options to normalize your data based on your study. You can also normalize your data using tools other than iCellR and import your data to iCellR. We recommend "ranked.glsf" normalization for most single cell studies. This normalization is great for fixing matrixes with lots of zeros and because it's geometric it will reduce some of batch differences in the library sizes, as long as all the data is aggregated into one file (to aggregate your data see "aggregating data" section above). GLSF stands for Geometric Library Size Factor, this is very similar to the normalization done by (DESeq2)[https://bioconductor.org/packages/release/bioc/html/DESeq2.html] and the ranked part would take the sum of the top most expressed genes as your library size instead of the full LB size which is to help resduce some of the drop out effects on normalization.
+## You have a few options to normalize your data based on your study. You can also normalize your data using tools other than iCellR and import your data to iCellR. We recommend "ranked.glsf" normalization for most single cell studies. This normalization is great for fixing matrixes with lots of zeros and because it's geometric it will reduce some of batch differences in the library sizes, as long as all the data is aggregated into one file (to aggregate your data see "aggregating data" section above). GLSF stands for Geometric Library Size Factor, this is very similar to the normalization done by [DESeq2](https://bioconductor.org/packages/release/bioc/html/DESeq2.html) and the ranked part would take the sum of the top most expressed genes as your library size instead of the full LB size which is to help resduce some of the drop out effects on normalization.
```r
my.obj <- norm.data(my.obj,
From a34f21cdf0bf0e59679016bbbd004fe07771ac1f Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:33:54 -0400
Subject: [PATCH 28/74] Update README.md
---
README.md | 4 +---
1 file changed, 1 insertion(+), 3 deletions(-)
diff --git a/README.md b/README.md
index f027fe1..b750532 100644
--- a/README.md
+++ b/README.md
@@ -306,9 +306,7 @@ my.obj <- norm.data(my.obj,
- Scale data (optional)
-iCellR dose not need this step as it scales the data when they need to be scaled on the fly; like for plotting or PCA.
-It is important to use the untansformed data for differential expression analysis to calculate the accurate fold changes.
-If you run this function the scaled data will be saved in different slot for you to download for plotting but will not be use by iCellR.
+iCellR does not need this step as it scales the data when they need to be scaled on the fly; like for plotting or running PCA. This is because, it is important to use the untransformed data for differential expression analysis to calculate the accurate/true fold changes. If you run this function the scaled data will replace the main data for this reason and instead will be saved in different data slot in the object for you to download if you need it for plotting or other reasons.
```r
# my.obj <- data.scale(my.obj)
From bc43a8b292fcd6ac15a9fd01d58519553e313f2e Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 11:35:59 -0400
Subject: [PATCH 29/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index b750532..6ecdf6a 100644
--- a/README.md
+++ b/README.md
@@ -336,7 +336,7 @@ head(my.obj@gene.data[order(my.obj@gene.data$numberOfCells, decreasing = T),])
- Make a gene model for clustering
-It's best to always to avoid global clustering and use a set of model genes. In bulk RNA-seq data it is very common to cluster the samples based on top 500 genes ranked by base mean, this is to reduce the noise. In scRNA-seq data, it's great to do so as well. This coupled with our ranked.glsf normalization is good for matrices with a lot of zeros. You can also use your set of genes as a model rather than making one.
+This function will help you find a good number of genes to use for running PCA.
```r
# See model plot
From 9cef3bc1e8525055dfde976c77781ac5b3d94090 Mon Sep 17 00:00:00 2001
From: Khodadadi-Jamayran
Date: Mon, 7 Jun 2021 12:01:48 -0400
Subject: [PATCH 30/74] 1.6.4
---
DESCRIPTION | 4 ++--
R/F00100.R | 2 +-
R/F0058.R | 21 +++++++++++++++------
R/F0071.R | 2 +-
man/i.score.Rd | 2 +-
5 files changed, 20 insertions(+), 11 deletions(-)
diff --git a/DESCRIPTION b/DESCRIPTION
index d544417..3a859e8 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -13,7 +13,7 @@ Authors@R: c(
person(given = "Aristotelis", family = "Tsirigos",role = c('aut','ctb'), comment = c(ORCID = '0000-0002-7512-8477'))
)
Maintainer: Alireza Khodadadi-Jamayran
-Description: A toolkit that allows scientists to work with data from single cell sequencing technologies such as scRNA-seq, scVDJ-seq and CITE-Seq. Single (i) Cell R package ('iCellR') provides unprecedented flexibility at every step of the analysis pipeline, including normalization, clustering, dimensionality reduction, imputation, visualization, and so on. Users can design both unsupervised and supervised models to best suit their research. In addition, the toolkit provides 2D and 3D interactive visualizations, differential expression analysis, filters based on cells, genes and clusters, data merging, normalizing for dropouts, data imputation methods, correcting for batch differences, pathway analysis, tools to find marker genes for clusters and conditions, predict cell types and pseudotime analysis. See Khodadadi-Jamayran, et al (2020) and Khodadadi-Jamayran, et al (2020) for more details.
+Description: A toolkit that allows scientists to work with data from single cell sequencing technologies such as scRNA-seq, scVDJ-seq, scATAC-seq, CITE-Seq and Spatial Transcriptomics (ST). Single (i) Cell R package ('iCellR') provides unprecedented flexibility at every step of the analysis pipeline, including normalization, clustering, dimensionality reduction, imputation, visualization, and so on. Users can design both unsupervised and supervised models to best suit their research. In addition, the toolkit provides 2D and 3D interactive visualizations, differential expression analysis, filters based on cells, genes and clusters, data merging, normalizing for dropouts, data imputation methods, correcting for batch differences, pathway analysis, tools to find marker genes for clusters and conditions, predict cell types and pseudotime analysis. See Khodadadi-Jamayran, et al (2020) and Khodadadi-Jamayran, et al (2020) for more details.
Depends: R (>= 3.3.0), ggplot2, plotly
Imports: Matrix, Rtsne, gridExtra, ggrepel, ggpubr, scatterplot3d,
RColorBrewer, knitr, NbClust, shiny, pheatmap, ape, ggdendro,
@@ -27,7 +27,7 @@ RoxygenNote: 7.1.1
BugReports: https://github.com/rezakj/iCellR/issues
URL: https://github.com/rezakj/iCellR
NeedsCompilation: yes
-Packaged: 2021-04-27 16:41:35 UTC; khodaa01
+Packaged: 2021-06-07 15:57:39 UTC; khodaa01
Author: Alireza Khodadadi-Jamayran [aut, cre]
(),
Joseph Pucella [aut, ctb] (),
diff --git a/R/F00100.R b/R/F00100.R
index 9c0ef8c..562ea29 100644
--- a/R/F00100.R
+++ b/R/F00100.R
@@ -162,7 +162,7 @@ utils::globalVariables(c("%>%",
# gg.cor.R | F0055.R
# globalVariables.R | F0056.R
# hto.anno.R | F0057.R
-# iNet.R | F0058.R
+# run.knetl.R | F0058.R
# iba.R | F0059.R
# iclust.R | F0060.R
# k.myImp.R | F0061.R
diff --git a/R/F0058.R b/R/F0058.R
index 19e8a68..48ee388 100644
--- a/R/F0058.R
+++ b/R/F0058.R
@@ -40,14 +40,22 @@ run.knetl <- function (x = NULL,
stop("x should be an object of class iCellR")
}
#####
- message("################# IMPORTANT NOTE ##################","")
+ message("############ IMPORTANT DISCLAIMER NOTE #############","")
message(paste(" Zoom is set to:", zoom,"! Make sure it is good for this data."))
message("
-# For data with less than 5000 cells use a zoom of about 100-200.
+ *** KNetL map is very dynamic with zoom and dims! ***
+ *** Therefore it needs to be adjusted! ***
+# For data with less than 1000 cells use a zoom of about 5-50.
+# For data with 1000-5000 cells use a zoom of about 50-200.
# For data with 5000-10000 cells use a zoom of about 100-300.
-# For data with 10000-30000 cells use a zoom of about 200-400.
+# For data with 10000-30000 cells use a zoom of about 200-500.
# For data with more than 30000 cells use a zoom of about 400-600.
-# A zoom value of 400 and dims = 1:20 is usually good for big data but adjust it for intended resolution.")
+# zoom 400 is usually good for big data but adjust for intended resolution.
+# Lower number for zoom in and higher for zoom out (its reverse).
+# dims = 1:20 is generally good for most data.
+# other parameters are best as default.
+
+ ")
message("###################################################","")
######
start_time1 <- Sys.time()
@@ -167,8 +175,6 @@ run.knetl <- function (x = NULL,
}
}
##################################
- # png('network_1.png',width = 30, height = 30, units = 'in', res = 300)
- #######
end_time1 <- Sys.time()
Time = difftime(end_time1,start_time1,units = "mins")
Time = round(as.numeric(Time),digits = 2)
@@ -217,3 +223,6 @@ run.knetl <- function (x = NULL,
return(g)
}
}
+
+
+
diff --git a/R/F0071.R b/R/F0071.R
index f6451a7..edb819f 100644
--- a/R/F0071.R
+++ b/R/F0071.R
@@ -5,7 +5,7 @@
#' @param data.type Choose from "raw.data" or "main.data", "imputed.data", default = "main.data".
#' @param scoring.List Genes that are used as a marker for phases.
#' @param return.stats Return the data or object. If FALSE the object would be returned.
-#' @param scoring.method Choose from "sum","mean" or "tirosh" for scoring method. See: https://science.sciencemag.org/content/352/6282/189
+#' @param scoring.method Choose from "tirosh (Tirosh, et. al. 2016), mean, sum, gsva, ssgsea, zscore and plage. , default = "tirosh". See: https://science.sciencemag.org/content/352/6282/189
#' @return The data frame object
#' @importFrom Hmisc cut2
#' @export
diff --git a/man/i.score.Rd b/man/i.score.Rd
index 8697f1f..f7ed5ea 100644
--- a/man/i.score.Rd
+++ b/man/i.score.Rd
@@ -21,7 +21,7 @@ i.score(
\item{return.stats}{Return the data or object. If FALSE the object would be returned.}
-\item{scoring.method}{Choose from "sum","mean" or "tirosh" for scoring method. See: https://science.sciencemag.org/content/352/6282/189}
+\item{scoring.method}{Choose from "tirosh (Tirosh, et. al. 2016), mean, sum, gsva, ssgsea, zscore and plage. , default = "tirosh". See: https://science.sciencemag.org/content/352/6282/189}
}
\value{
The data frame object
From fad2862933cf335ca053244ae32c66faa2742c1c Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 8 Jun 2021 11:58:57 -0400
Subject: [PATCH 31/74] Update .travis.yml
---
.travis.yml | 25 +++++++++++++++++++++++++
1 file changed, 25 insertions(+)
diff --git a/.travis.yml b/.travis.yml
index 8d240c7..e363703 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -5,6 +5,7 @@ sudo: required
cache: packages
os:
+ - linux
- osx
r:
@@ -16,4 +17,28 @@ matrix:
- r: devel
os: osx
+env:
+ global:
+ - _R_CHECK_FORCE_SUGGESTS_=FALSE
+ - ASAN="-fsanitize=address -fno-omit-frame-pointer"
+ - LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/usr/local/lib
+ - HDF5_VERSION=1.8.17
+ - HDF5_RELEASE_URL="https://support.hdfgroup.org/ftp/HDF5/releases"
+
+before_install:
+ - chmod +x travis_setup.sh
+ - ./travis_setup.sh
+
+addons:
+ apt:
+ packages:
+ - subversion
+ - autoconf
+ - build-essential
+ - libtool
+ - libmagick++-dev
+ homebrew:
+ packages:
+ - libgit2
+
warnings_are_errors: false
From 85a8a034e4e21b57615fd59a3f8cf0522f424817 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 8 Jun 2021 12:02:41 -0400
Subject: [PATCH 32/74] Delete .travis.yml
---
.travis.yml | 44 --------------------------------------------
1 file changed, 44 deletions(-)
delete mode 100644 .travis.yml
diff --git a/.travis.yml b/.travis.yml
deleted file mode 100644
index e363703..0000000
--- a/.travis.yml
+++ /dev/null
@@ -1,44 +0,0 @@
-# R for travis: see documentation at https://docs.travis-ci.com/user/languages/r
-
-language: r
-sudo: required
-cache: packages
-
-os:
- - linux
- - osx
-
-r:
- - release
- - devel
-
-matrix:
- exclude:
- - r: devel
- os: osx
-
-env:
- global:
- - _R_CHECK_FORCE_SUGGESTS_=FALSE
- - ASAN="-fsanitize=address -fno-omit-frame-pointer"
- - LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/usr/local/lib
- - HDF5_VERSION=1.8.17
- - HDF5_RELEASE_URL="https://support.hdfgroup.org/ftp/HDF5/releases"
-
-before_install:
- - chmod +x travis_setup.sh
- - ./travis_setup.sh
-
-addons:
- apt:
- packages:
- - subversion
- - autoconf
- - build-essential
- - libtool
- - libmagick++-dev
- homebrew:
- packages:
- - libgit2
-
-warnings_are_errors: false
From 89cb4ec06fc9adbb7b7e01a297fbc121082528ef Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 8 Jun 2021 12:02:53 -0400
Subject: [PATCH 33/74] Update README.md
---
README.md | 1 -
1 file changed, 1 deletion(-)
diff --git a/README.md b/README.md
index 6ecdf6a..4b98c2b 100644
--- a/README.md
+++ b/README.md
@@ -1,7 +1,6 @@
[![CRAN Version](https://www.r-pkg.org/badges/version/iCellR)](https://cran.r-project.org/package=iCellR)
[![CRAN Downloads](https://cranlogs.r-pkg.org/badges/iCellR)](https://cran.r-project.org/package=iCellR)
[![License: GPL v2](https://img.shields.io/badge/License-GPL%20v2-blue.svg)](https://www.gnu.org/licenses/old-licenses/gpl-2.0.en.html)
-[![Build Status](https://travis-ci.com/rezakj/iCellR.svg?branch=master)](https://travis-ci.com/rezakj/iCellR)
# iCellR
iCellR is an interactive R package to work with high-throughput single cell sequencing technologies (i.e scRNA-seq, scVDJ-seq, scATAC-seq, CITE-Seq and Spatial Transcriptomics (ST)).
From fcbdeb8802c9cc7723f070bc421c9049e452820e Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Wed, 9 Jun 2021 12:19:00 -0400
Subject: [PATCH 34/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 4b98c2b..d5dfcea 100644
--- a/README.md
+++ b/README.md
@@ -505,7 +505,7 @@ In this function we provide a variety of many other options for you to explore t
| ------------- | ------------- | ------------- |
| ward.D, ward.D2, single, complete, average, mcquitty, median, centroid, kmeans| euclidean, maximum, manhattan, canberra, binary, minkowski or NULL | kl, ch, hartigan, ccc, scott, marriot, trcovw, tracew, friedman, rubin, cindex, db, silhouette, duda, pseudot2, beale, ratkowsky, ball, ptbiserial, gap, frey, mcclain, gamma, gplus, tau, dunn, hubert, sdindex, dindex, sdbw |
-### Conventionally people cluster based on PCA data (usually first 10 dimensions) however you have the option of choosing tSNA, UMAP and KNetL map as well. If you have adjusted your KNetL map and are confident about the results we recommend clustering based on KNetL map.
+### Conventionally people cluster based on PCA data (usually first 10 dimensions) however you have the option of choosing tSNE, UMAP and KNetL map dimensions as well. If you have adjusted your KNetL map and are confident about the results we recommend clustering based on KNetL map.
This is one of the harder parts of the analysis and sometimes you need to adjust your clustering based on marker genes. This means you might need to merge some clusters, gate (see our cell gating tools) or try different sensitivities to find more or less communities.
From 66f733756f1f9272d83d3e282c02df9bc4755c70 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Mon, 9 Aug 2021 16:44:20 -0400
Subject: [PATCH 35/74] Update README.md
---
README.md | 7 +++----
1 file changed, 3 insertions(+), 4 deletions(-)
diff --git a/README.md b/README.md
index d5dfcea..066272c 100644
--- a/README.md
+++ b/README.md
@@ -23,10 +23,9 @@ iCellR Viewer (web GUI app): https://compbio.nyumc.org/icellr/
If you are using FlowJo or SeqGeq, they have made plugins for iCellR and other single cell tools: https://www.flowjo.com/exchange/#/ (list of all plugins) and https://www.flowjo.com/exchange/#/plugin/profile?id=34 (iCellR plugin). [SeqGeq DE tutorial](https://www.youtube.com/watch?v=gXFmWRpdwow)
-For citing iCellR use these:
-- Citation for batch alignment and imputation: https://www.biorxiv.org/content/10.1101/2020.03.31.019109v1.full
-- Citation for KNetL map: https://www.biorxiv.org/content/10.1101/2020.05.05.078550v1.full
-- iCellR publications: [PMID 31744829](https://www.ncbi.nlm.nih.gov/pubmed/31744829) (scRNA-seq), [PMID: 31934613](https://www.ncbi.nlm.nih.gov/pubmed/31934613) (bulk RNA-seq from TCGA), [PMID: 32550269](https://pubmed.ncbi.nlm.nih.gov/32550269/) (scVDJ-seq)
+For citing iCellR use this [PMID: 34353854](https://cancerdiscovery.aacrjournals.org/content/early/2021/07/28/2159-8290.CD-21-0369)
+
+iCellR publications: [PMID 31744829](https://www.ncbi.nlm.nih.gov/pubmed/31744829) (scRNA-seq), [PMID: 31934613](https://www.ncbi.nlm.nih.gov/pubmed/31934613) (bulk RNA-seq from TCGA), [PMID: 32550269](https://pubmed.ncbi.nlm.nih.gov/32550269/) (scVDJ-seq), [PMID: 34135081](https://jasn.asnjournals.org/content/32/8/1987), [PMID: 33593073](https://www.ahajournals.org/doi/10.1161/CIRCRESAHA.120.317914)
### Single (i) Cell R package (iCellR)
From 879ac988504c27c13d2ea0f33d76022f5db90989 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Thu, 26 Aug 2021 15:42:02 -0400
Subject: [PATCH 36/74] Update README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 066272c..e4ec8e1 100644
--- a/README.md
+++ b/README.md
@@ -5,7 +5,7 @@
# iCellR
iCellR is an interactive R package to work with high-throughput single cell sequencing technologies (i.e scRNA-seq, scVDJ-seq, scATAC-seq, CITE-Seq and Spatial Transcriptomics (ST)).
-### News (April 2021): Use iCellR version 1.6.4 for scATAC-seq and Spatial Transcriptomics (ST). Use function i.score for scoring (scoring cells based on genes of your interest) methods (i.e. tirosh, mean, sum, gsva, ssgsea, zscore and plage).
+### News (April 2021): Use iCellR version 1.6.4 for scATAC-seq and Spatial Transcriptomics (ST). Use the i.score function for scoring (scoring cells based on gene signatures) methods (i.e. tirosh, mean, sum, gsva, ssgsea, zscore and plage).
### News (July 2020): See iCellR version 1.5.5 with new cell cycle analysis for G0, G1S, G2M, M, G1M and S [phase](https://genome.med.nyu.edu/results/external/iCellR/example1/All_cellcycle.png), Pseudotime Abstract KNetL map [(PAK map)](https://genome.med.nyu.edu/results/external/iCellR/example1/pseudotime.KNetL.png) and gene-gene [correlations](https://genome.med.nyu.edu/results/external/iCellR/example1/gene-gene.correlation.png). See below for how to.
From 49cc5c43c78b9312a490c2bcc1b714366645c98f Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 7 Sep 2021 18:33:47 -0400
Subject: [PATCH 37/74] Update README.md
---
README.md | 44 ++++++++++++++++++++++++++++++++++++++++++--
1 file changed, 42 insertions(+), 2 deletions(-)
diff --git a/README.md b/README.md
index e4ec8e1..987a37d 100644
--- a/README.md
+++ b/README.md
@@ -2931,7 +2931,7 @@ heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters"
dev.off()
```
-Work on scATAC data
+Work on scATAC data (normalize and find marker peaks for each cluster)
```r
# normalize ACAT
@@ -3004,12 +3004,34 @@ my.obj <- run.impute(my.obj,data.type = "knetl", nn = 10, ATAC.data = FALSE)
png('heatmap_gg_peaks.png', width = 10, height = 10, units = 'in', res = 300)
heatmap.gg.plot(my.obj, gene = MyGenes, interactive = F, cluster.by = "clusters",cell.sort = F, conds.to.plot = NULL, data.type = "atac.imputed")
dev.off()
+
+## you can also find avarage peak intensity per cluster
+
+my.obj <- clust.avg.exp(my.obj, data.type = "atac")
+head(my.obj@clust.avg)
+
+#gene cluster_1 cluster_2 ...
+#chr1.100037799.100038931 0.38238731 0.36750000 ...
+#chr1.100132733.100133298 0.11195725 1.13593827 ...
+#chr1.100249637.100250160 0.09851425 0.09511728 ...
+#chr1.100265992.100266479 0.06768394 0.17707407 ...
+#chr1.10032488.10033387 0.35273705 0.14885802 ...
+#chr1.100352150.100352921 0.12006088 0.00000000 ...
+
+# find out which cluster has the highest number
+
+dat <- as.data.frame(t((my.obj@clust.avg)[,-1]))
+dat <- hto.anno(hto.data = dat)
+
+head(dat$assignment.annotatio)
+#[1] cluster_1 cluster_2 cluster_4 cluster_2 cluster_3 cluster_4
+#8 Levels: cluster_1 cluster_2 cluster_3 cluster_4 cluster_5 ... cluster_8
```
Peak analysis
```r
-# make bed file per cluster
+# make a bed file per cluster from the marker.peaks file you made up here
make.bed(marker.peaks)
# load packages
@@ -3070,8 +3092,26 @@ lapply(1:length(genes), function(i) write.table(genes[[i]],
```
+Merging scATAC files with different intervals (as dipicted in bedtools website)
+
+
+
+
+
+```
+# Let's say you have 3 files that you need to merege
+# example file
+head(File1)[1:3]
+# AAACAGCCAAGTGAAC.1 AAACAGCCACTGACCG.1 AAACAGCCATGATTGT.1
+#chr1.181218.181695 0 0 1
+#chr1.191296.191699 0 0 0
+#chr1.629770.630129 0 0 0
+#chr1.633806.634251 0 0 0
+#chr1.778422.779040 0 0 0
+#chr1.827306.827702 0 0 0
+```
From 3d4e206febe16b365df54621fe14d191cb34f771 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 7 Sep 2021 19:20:42 -0400
Subject: [PATCH 38/74] Update README.md
---
README.md | 131 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
1 file changed, 131 insertions(+)
diff --git a/README.md b/README.md
index 987a37d..be6ac19 100644
--- a/README.md
+++ b/README.md
@@ -3111,6 +3111,137 @@ head(File1)[1:3]
#chr1.778422.779040 0 0 0
#chr1.827306.827702 0 0 0
+# get the row names from each file and concatenate them as below:
+
+f1 <- row.names(File1)
+f2 <- row.names(File2)
+f3 <- row.names(File3)
+
+all.peaks <- c(f1,f2,f3)
+head(all.peaks)
+#[1] "chr1.181218.181695" "chr1.191296.191699" "chr1.629770.630129"
+#[4] "chr1.633806.634251" "chr1.778422.779040" "chr1.827306.827702"
+
+chr <- as.character(as.matrix(data.frame(do.call('rbind', strsplit(as.character(all.peaks),'.',fixed=TRUE)))[1]))
+start <- data.frame(do.call('rbind', strsplit(as.character(all.peaks),'.',fixed=TRUE)))[2]
+end <- data.frame(do.call('rbind', strsplit(as.character(all.peaks),'.',fixed=TRUE)))[3]
+
+# make a bed file
+
+DAT <- as.data.frame(chr)
+DAT$start <- as.numeric(as.matrix(start))
+DAT$end <- as.numeric(as.matrix(end))
+head(DAT)
+# chr start end
+#1 chr1 181218 181695
+#2 chr1 191296 191699
+#3 chr1 629770 630129
+#4 chr1 633806 634251
+#5 chr1 778422 779040
+#6 chr1 827306 827702
+
+# make a GenomicRanges object
+
+library("GenomicRanges")
+
+all.gr <- GRanges(seqnames=DAT$chr,ranges=IRanges(start=DAT$start,end=DAT$end))
+
+all.gr
+#GRanges object with ?? ranges and 0 metadata columns:
+# seqnames ranges strand
+#
+# [1] chr1 181218-181695 *
+# [2] chr1 191296-191699 *
+# [3] chr1 629770-630129 *
+# [4] chr1 633806-634251 *
+# [5] chr1 778422-779040 *
+# ... ... ... ...
+# [52] chr1 1303892-1306216 *
+# [53] chr1 1307242-1309359 *
+# [54] chr1 1324425-1325236 *
+# [55] chr1 1348940-1349958 *
+# [56] chr1 1372031-1372220 *
+# -------
+# seqinfo: 1 sequence from an unspecified genome; no seqlengths
+
+# sort and merge the peaks
+
+mrg <- reduce(all.gr)
+
+########################## choose file and give name
+
+MyFile <- f1
+name="f1_new.bed"
+
+########################## cop paste the code here to make a new bed file that has the old and new intervals to be replaced
+
+chr <- as.character(as.matrix(data.frame(do.call('rbind', strsplit(as.character(MyFile),'.',fixed=TRUE)))[1]))
+start <- data.frame(do.call('rbind', strsplit(as.character(MyFile),'.',fixed=TRUE)))[2]
+end <- data.frame(do.call('rbind', strsplit(as.character(MyFile),'.',fixed=TRUE)))[3]
+# make a bed file
+DAT <- as.data.frame(chr)
+DAT$start <- as.numeric(as.matrix(start))
+DAT$end <- as.numeric(as.matrix(end))
+MyFile <- DAT
+
+# make intrval file to replace to new regions
+MyFile.gr <- GRanges(seqnames=MyFile$chr,ranges=IRanges(start=MyFile$start,end=MyFile$end))
+
+F <- subsetByOverlaps(MyFile.gr,mrg)
+M <- subsetByOverlaps(mrg,MyFile.gr)
+###
+chr <- as.character(F@seqnames)
+DAT <- as.data.frame(chr)
+DAT$start <- F@ranges@start
+DAT$end <- (F@ranges@start + F@ranges@width) - 1
+DAT$new.chr<- as.character(M@seqnames)
+DAT$new.start <- M@ranges@start
+DAT$new.end <- (M@ranges@start + M@ranges@width) - 1
+# diff
+ADD <- setdiff(mrg,MyFile.gr)
+L <- length(as.character(ADD@seqnames))
+chr <-rep("NA",L)
+DAT1 <- as.data.frame(chr)
+DAT1$start <- rep("NA",L)
+DAT1$end <- rep("NA",L)
+DAT1$new.chr<- as.character(ADD@seqnames)
+DAT1$new.start <- ADD@ranges@start
+DAT1$new.end <- (ADD@ranges@start + ADD@ranges@width) - 1
+
+Final.DAT <- rbind(DAT,DAT1)
+
+##### Write
+
+write.table(Final.DAT,name,row.names=FALSE,sep="\t", quote = FALSE)
+
+# example file
+# head(Final.DAT,10)
+# chr start end new.chr new.start new.end
+#1 chr1 181218 181695 chr1 181218 181695
+#2 chr1 191296 191699 chr1 191296 191699
+#3 chr1 629770 630129 chr1 629770 630129
+#4 chr1 633806 634251 chr1 633806 634251
+#5 chr1 778422 779040 chr1 778422 779040
+#6 chr1 827306 827702 chr1 827306 827702
+#7 NA NA NA chr1 904635 904943
+#8 NA NA NA chr1 923684 924085
+#9 NA NA NA chr1 940590 940783
+#10 NA NA NA chr1 959052 959594
+
+# The first 3 columns are the original peaks and the last 3 are the ones that need to be replaced with original one. The NA peaks would also get the new peak ids but in the matrix the cells will have 0 for expressions. To do this use the iCellR function replace.peak.id.
+
+MyATAC1 <- replace.peak.id(atac.data=MyATAC1, bed.file = Final.DAT)
+MyATAC2 <- replace.peak.id(atac.data=MyATAC2, bed.file = Final.DAT)
+MyATAC3 <- replace.peak.id(atac.data=MyATAC3, bed.file = Final.DAT)
+
+# finally aggregate the samples and add to iCellR object
+
+my.atac.data <- data.aggregation(samples = c("MyATAC1","MyATAC2","MyATAC3"),
+ condition.names = c("WT","KO","Ctrl"))
+
+# add ATAC-Seq data
+my.obj@atac.raw <- MyATAC
+my.obj@atac.main <- MyATAC
```
From 884dfee2b3ce08e6f4679ad6e9abc93b864048a2 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Tue, 7 Sep 2021 19:26:27 -0400
Subject: [PATCH 39/74] Update README.md
---
README.md | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/README.md b/README.md
index be6ac19..21e7b65 100644
--- a/README.md
+++ b/README.md
@@ -3240,8 +3240,8 @@ my.atac.data <- data.aggregation(samples = c("MyATAC1","MyATAC2","MyATAC3"),
condition.names = c("WT","KO","Ctrl"))
# add ATAC-Seq data
-my.obj@atac.raw <- MyATAC
-my.obj@atac.main <- MyATAC
+my.obj@atac.raw <- my.atac.data
+my.obj@atac.main <- my.atac.data
```
From 19a825469735003fc783a5018d7466d65aecf617 Mon Sep 17 00:00:00 2001
From: Alireza Khodadadi-Jamayran
Date: Thu, 9 Sep 2021 11:09:41 -0400
Subject: [PATCH 40/74] Update README.md
---
README.md | 3 ++-
1 file changed, 2 insertions(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 21e7b65..08d0caf 100644
--- a/README.md
+++ b/README.md
@@ -3173,7 +3173,8 @@ mrg <- reduce(all.gr)
MyFile <- f1
name="f1_new.bed"
-########################## cop paste the code here to make a new bed file that has the old and new intervals to be replaced
+########################## copy paste the code here to make a new bed file
+########################## the new bed has the old and new intervals (new intervals to be replaced with old)
chr <- as.character(as.matrix(data.frame(do.call('rbind', strsplit(as.character(MyFile),'.',fixed=TRUE)))[1]))
start <- data.frame(do.call('rbind', strsplit(as.character(MyFile),'.',fixed=TRUE)))[2]
From 71ca6785eddd74968e1c83006b89ba617d2086aa Mon Sep 17 00:00:00 2001
From: Khodadadi-Jamayran
Date: Thu, 7 Oct 2021 14:33:30 -0400
Subject: [PATCH 41/74] =?UTF-8?q?=E2=80=9C1.6.5=E2=80=9D?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
---
DESCRIPTION | 4 ++--
R/{F001.R => F0001.R} | 0
R/{F002.R => F0002.R} | 0
R/{F003.R => F0003.R} | 0
R/{F004.R => F0004.R} | 3 ---
R/{F005.R => F0005.R} | 10 ----------
R/{F006.R => F0006.R} | 11 -----------
R/{F007.R => F0007.R} | 5 -----
R/{F008.R => F0008.R} | 4 ----
R/{F009.R => F0009.R} | 5 -----
R/F0010.R | 4 ----
R/F0011.R | 20 --------------------
R/F0012.R | 5 -----
R/F0013.R | 3 ---
R/F0014.R | 6 ------
R/F0015.R | 11 -----------
R/F0017.R | 5 -----
R/F0018.R | 4 ----
R/F0020.R | 12 ------------
R/F0025.R | 9 ---------
R/F0027.R | 5 -----
R/F0029.R | 3 ---
R/F0030.R | 10 ----------
R/F0031.R | 15 ---------------
R/F0032.R | 10 ----------
R/F0034.R | 5 -----
R/F0035.R | 16 ----------------
R/F0039.R | 11 -----------
R/F0040.R | 13 -------------
R/F0042.R | 13 -------------
R/F0043.R | 5 -----
R/F0044.R | 1 +
R/F0052.R | 31 +++++++++++++++++++++++++------
R/F0057.R | 11 -----------
R/F0072.R | 26 --------------------------
R/{F00100.R => F0100.R} | 4 ++--
data/demo.obj.rda | Bin 257086 -> 0 bytes
man/add.vdj.Rd | 14 --------------
man/cell.filter.Rd | 14 +-------------
man/change.clust.Rd | 10 ----------
man/clust.stats.plot.Rd | 6 ------
man/cluster.plot.Rd | 13 -------------
man/data.aggregation.Rd | 2 +-
man/data.scale.Rd | 8 +-------
man/demo.obj.Rd | 19 -------------------
man/down.sample.Rd | 8 +-------
man/find.dim.genes.Rd | 7 -------
man/gene.plot.Rd | 16 ----------------
man/gene.stats.Rd | 5 -----
man/heatmap.gg.plot.Rd | 11 -----------
man/hto.anno.Rd | 12 ------------
man/load10x.Rd | 2 +-
man/make.bed.Rd | 2 +-
man/make.gene.model.Rd | 21 ---------------------
man/make.obj.Rd | 2 +-
man/norm.data.Rd | 7 +------
man/opt.pcs.plot.Rd | 4 ----
man/prep.vdj.Rd | 12 ------------
man/pseudotime.tree.Rd | 11 -----------
man/qc.stats.Rd | 6 +-----
man/run.clustering.Rd | 12 ------------
man/run.diff.exp.Rd | 6 ------
man/run.pc.tsne.Rd | 6 ------
man/run.pca.Rd | 6 ------
man/run.tsne.Rd | 6 ------
man/run.umap.Rd | 5 -----
man/stats.plot.Rd | 13 +------------
man/top.markers.Rd | 4 ----
man/vdj.stats.Rd | 14 --------------
man/volcano.ma.plot.Rd | 17 -----------------
70 files changed, 40 insertions(+), 561 deletions(-)
rename R/{F001.R => F0001.R} (100%)
rename R/{F002.R => F0002.R} (100%)
rename R/{F003.R => F0003.R} (100%)
rename R/{F004.R => F0004.R} (97%)
rename R/{F005.R => F0005.R} (96%)
rename R/{F006.R => F0006.R} (95%)
rename R/{F007.R => F0007.R} (97%)
rename R/{F008.R => F0008.R} (97%)
rename R/{F009.R => F0009.R} (85%)
delete mode 100644 R/F0072.R
rename R/{F00100.R => F0100.R} (98%)
delete mode 100644 data/demo.obj.rda
delete mode 100644 man/demo.obj.Rd
diff --git a/DESCRIPTION b/DESCRIPTION
index 3a859e8..5c22100 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
Package: iCellR
Type: Package
Title: Analyzing High-Throughput Single Cell Sequencing Data
-Version: 1.6.4
+Version: 1.6.5
Authors@R: c(
person(given = 'Alireza', family = 'Khodadadi-Jamayran',role = c('aut','cre'), email = 'alireza.khodadadi.j@gmail.com', comment = c(ORCID = '0000-0003-2495-7504')),
person(given = 'Joseph', family = 'Pucella',role = c('aut','ctb'), comment = c(ORCID = '0000-0003-0875-8046')),
@@ -27,7 +27,7 @@ RoxygenNote: 7.1.1
BugReports: https://github.com/rezakj/iCellR/issues
URL: https://github.com/rezakj/iCellR
NeedsCompilation: yes
-Packaged: 2021-06-07 15:57:39 UTC; khodaa01
+Packaged: 2021-10-07 18:26:50 UTC; khodaa01
Author: Alireza Khodadadi-Jamayran [aut, cre]
(),
Joseph Pucella [aut, ctb] (),
diff --git a/R/F001.R b/R/F0001.R
similarity index 100%
rename from R/F001.R
rename to R/F0001.R
diff --git a/R/F002.R b/R/F0002.R
similarity index 100%
rename from R/F002.R
rename to R/F0002.R
diff --git a/R/F003.R b/R/F0003.R
similarity index 100%
rename from R/F003.R
rename to R/F0003.R
diff --git a/R/F004.R b/R/F0004.R
similarity index 97%
rename from R/F004.R
rename to R/F0004.R
index 4c0e37e..2f0e368 100644
--- a/R/F004.R
+++ b/R/F0004.R
@@ -7,9 +7,6 @@
#' @param s.phase.genes A character vector of gene names for S phase, default = s.phase.
#' @param g2m.phase.genes A character vector of gene names for G2 and M phase, default = g2m.phase.
#' @return The data frame object
-#' @examples
-#' New.demo.obj <- qc.stats(demo.obj)
-#' head(New.demo.obj@stats)
#' @export
qc.stats <- function (x = NULL,
which.data = "raw.data",
diff --git a/R/F005.R b/R/F0005.R
similarity index 96%
rename from R/F005.R
rename to R/F0005.R
index 06ebd35..17bd4d2 100644
--- a/R/F005.R
+++ b/R/F0005.R
@@ -12,16 +12,6 @@
#' @param interactive If set to TRUE an interactive HTML file will be created, default = TRUE.
#' @param out.name If "interactive" is set to TRUE, the out put name for HTML, default = "plot".
#' @return An object of class iCellR.
-#' @examples
-#' stats.plot(demo.obj,
-#' plot.type = "three.in.one",
-#' out.name = "UMI-plot",
-#' interactive = FALSE,
-#' cell.color = "slategray3",
-#' cell.size = 1,
-#' cell.transparency = 0.5,
-#' box.color = "red",
-#' box.line.col = "green")
#' @import gridExtra
#' @importFrom htmlwidgets saveWidget
#' @importFrom plotly ggplotly layout plot_ly
diff --git a/R/F006.R b/R/F0006.R
similarity index 95%
rename from R/F006.R
rename to R/F0006.R
index cd43021..26de36c 100644
--- a/R/F006.R
+++ b/R/F0006.R
@@ -13,17 +13,6 @@
#' @param filter.by.gene A character vector of gene names to be filtered by thier expression. If more then one gene is defined it would be OR not AND.
#' @param filter.by.gene.exp.min Minimum gene expression to be filtered by the genes set in filter.by.gene, default = 1.
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- cell.filter(demo.obj,
-#' min.mito = 0,
-#' max.mito = 0.05 ,
-#' min.genes = 100,
-#' max.genes = 2500,
-#' min.umis = 0,
-#' max.umis = Inf)
-#'
-#' message(demo.obj@my.filters)
-#'
#' @export
cell.filter <- function (x = NULL,
min.mito = 0,
diff --git a/R/F007.R b/R/F0007.R
similarity index 97%
rename from R/F007.R
rename to R/F0007.R
index d6a35a7..8428192 100644
--- a/R/F007.R
+++ b/R/F0007.R
@@ -3,11 +3,6 @@
#' This function takes an object of class iCellR and down samples the condition to have equal number of cells in each condition.
#' @param x An object of class iCellR.
#' @return An object of class iCellR.
-#' @examples
-#'
-#' my.obj <- down.sample(demo.obj)
-#'
-#'
#' @export
down.sample <- function (x = NULL) {
if ("iCellR" != class(x)[1]) {
diff --git a/R/F008.R b/R/F0008.R
similarity index 97%
rename from R/F008.R
rename to R/F0008.R
index c5e2125..8f7474a 100644
--- a/R/F008.R
+++ b/R/F0008.R
@@ -12,10 +12,6 @@
#' @param ATAC.data If TURE, it would normalize ATAC-Seq data and not RNA-Seq, default = FALSE.
#' @param ATAC.filter If TURE, all the cells filtered in RNA-Seq will be filtered in ATAC-Seq. This needs to be done for both data to match, default = TRUE.
#' @return An object of class iCellR.
-#' @examples
-#'
-#' demo.obj <- norm.data(demo.obj, norm.method = "ranked.glsf", top.rank = 500)
-#'
#' @export
norm.data <- function (x = NULL,
norm.method = "ranked.glsf",
diff --git a/R/F009.R b/R/F0009.R
similarity index 85%
rename from R/F009.R
rename to R/F0009.R
index 33422cf..91403d2 100644
--- a/R/F009.R
+++ b/R/F0009.R
@@ -3,11 +3,6 @@
#' This function takes an object of class iCellR and scales the normalized data.
#' @param x An object of class iCellR.
#' @return An object of class iCellR.
-#' @examples
-#' my.obj <- data.scale(demo.obj)
-#'
-#' head(my.obj@scaled.data)[1:5]
-#'
#' @export
data.scale <- function (x = NULL) {
if ("iCellR" != class(x)[1]) {
diff --git a/R/F0010.R b/R/F0010.R
index 0ea8172..4b08081 100644
--- a/R/F0010.R
+++ b/R/F0010.R
@@ -5,10 +5,6 @@
#' @param which.data Choose from "raw.data" or "main.data", default = "raw.data".
#' @param each.cond If TRUE each condition will be calculated, default = FALSE.
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- gene.stats(demo.obj, which.data = "main.data")
-#' head(demo.obj@gene.data)
-#'
#' @export
gene.stats <- function (x = NULL,
which.data = "raw.data",
diff --git a/R/F0011.R b/R/F0011.R
index f82d782..4f5b5c0 100644
--- a/R/F0011.R
+++ b/R/F0011.R
@@ -19,26 +19,6 @@
#' @param interactive If set to TRUE an interactive HTML file will be created, default = TRUE.
#' @param out.name If "interactive" is set to TRUE, the out put name for HTML, default = "plot".
#' @return An object of class iCellR.
-#' @examples
-#' make.gene.model(demo.obj,
-#' dispersion.limit = 1.5,
-#' base.mean.rank = 500,
-#' no.mito.model = TRUE,
-#' mark.mito = TRUE,
-#' interactive = FALSE,
-#' my.out.put = "plot",
-#' out.name = "gene.model")
-#'
-#' demo.obj <- make.gene.model(demo.obj,
-#' dispersion.limit = 1.5,
-#' base.mean.rank = 500,
-#' no.mito.model = TRUE,
-#' mark.mito = TRUE,
-#' interactive = FALSE,
-#' out.name = "gene.model")
-#'
-#' head(demo.obj@gene.model)
-#'
#' @import ggrepel
#' @export
make.gene.model <- function (x = NULL,
diff --git a/R/F0012.R b/R/F0012.R
index 30636cd..a9adaf5 100644
--- a/R/F0012.R
+++ b/R/F0012.R
@@ -9,11 +9,6 @@
#' @param gene.list A charactor vector of genes to be used for PCA. If "clust.method" is set to "gene.model", default = "my_model_genes.txt".
#' @param scale.data If TRUE the data will be scaled (log2 + plus.log.value), default = TRUE.
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- run.pca(demo.obj, method = "gene.model", gene.list = demo.obj@gene.model)
-#'
-#' head(demo.obj@pca.data)[1:5]
-#'
#' @export
run.pca <- function (x = NULL,
data.type = "main",
diff --git a/R/F0013.R b/R/F0013.R
index 9780a65..583be6e 100644
--- a/R/F0013.R
+++ b/R/F0013.R
@@ -4,9 +4,6 @@
#' @param x An object of class iCellR.
#' @param pcs.in.plot Number of PCs to show in plot, defult = 50.
#' @return An object of class iCellR.
-#' @examples
-#' opt.pcs.plot(demo.obj)
-#'
#' @export
opt.pcs.plot <- function (x = NULL, pcs.in.plot = 50) {
if ("iCellR" != class(x)[1]) {
diff --git a/R/F0014.R b/R/F0014.R
index e1a734f..f72262f 100644
--- a/R/F0014.R
+++ b/R/F0014.R
@@ -6,12 +6,6 @@
#' @param top.pos Number of top positive marker genes to be taken from each PC, default = 15.
#' @param top.neg Number of top negative marker genes to be taken from each PC, default = 5.
#' @return An object of class iCellR.
-#' @examples
-#'
-#' demo.obj <- find.dim.genes(demo.obj, dims = 1:10,top.pos = 20, top.neg = 20)
-#'
-#' head(demo.obj@gene.model)
-#'
#' @export
find.dim.genes <- function (x = NULL,
dims = 1:10,
diff --git a/R/F0015.R b/R/F0015.R
index 482c73b..8e21eba 100644
--- a/R/F0015.R
+++ b/R/F0015.R
@@ -9,17 +9,6 @@
#' @param min.clust minimum number of clusters, default = 2.
#' @param dims PCA dimentions to be use for clustering, default = 1:10.
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- run.clustering(demo.obj,
-#' clust.method = "kmeans",
-#' dist.method = "euclidean",
-#' index.method = "silhouette",
-#' max.clust = 2,
-#' min.clust = 2,
-#' dims = 1:10)
-#'
-#' head(demo.obj@best.clust)
-#'
#' @import NbClust
#' @export
run.clustering <- function (x = NULL,
diff --git a/R/F0017.R b/R/F0017.R
index eb6136e..6047edf 100644
--- a/R/F0017.R
+++ b/R/F0017.R
@@ -23,11 +23,6 @@
#' @param eta numeric; Learning rate (default: 200.0)
#' @param exaggeration_factor numeric; Exaggeration factor used to multiply the P matrix in the first part of the optimization (default: 12.0)
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- run.pc.tsne(demo.obj, dims = 1:10,perplexity = 20)
-#'
-#' head(demo.obj@pca.data)[1:5]
-#'
#' @import Rtsne
#' @export
run.pc.tsne <- function (x = NULL,
diff --git a/R/F0018.R b/R/F0018.R
index 015a9f6..342d8f1 100644
--- a/R/F0018.R
+++ b/R/F0018.R
@@ -286,10 +286,6 @@
#' @param verbose If \code{TRUE}, log details to the console.
#'
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- run.umap(demo.obj, dims = 1:10)
-#' head(demo.obj@umap.data)
-#'
#' @import uwot
#' @export
run.umap <- function (x = NULL,
diff --git a/R/F0020.R b/R/F0020.R
index 76064ac..c50e780 100644
--- a/R/F0020.R
+++ b/R/F0020.R
@@ -20,18 +20,6 @@
#' @param density If TRUE the density plots for PCA/tSNE second dimension will be created, default = FALSE.
#' @param static3D If TRUE a non-interactive 3D plot will be made.
#' @return An object of class iCellR.
-#' @examples
-#' cluster.plot(demo.obj,plot.type = "umap",interactive = FALSE)
-#'
-#' cluster.plot(demo.obj,plot.type = "tsne",interactive = FALSE)
-#'
-#' cluster.plot(demo.obj,plot.type = "pca",interactive = FALSE)
-#'
-#' cluster.plot(demo.obj,plot.type = "pca",col.by = "conditions",interactive = FALSE)
-#'
-#' cluster.plot(demo.obj,plot.type = "umap",col.by = "conditions",interactive = FALSE)
-#'
-#' cluster.plot(demo.obj,plot.type = "tsne",col.by = "conditions",interactive = FALSE)
#' @import RColorBrewer
#' @import scatterplot3d
#' @importFrom htmlwidgets saveWidget
diff --git a/R/F0025.R b/R/F0025.R
index 9697331..758e169 100644
--- a/R/F0025.R
+++ b/R/F0025.R
@@ -6,15 +6,6 @@
#' @param to.clust The new name for the cluster.
#' @param clust.reset Reset to the original clustering.
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- change.clust(demo.obj, change.clust = 1, to.clust = 3)
-#' cluster.plot(demo.obj,plot.type = "umap",interactive = FALSE)
-#'
-#' demo.obj <- change.clust(demo.obj, change.clust = 3, to.clust = "B Cell")
-#' cluster.plot(demo.obj,plot.type = "umap",interactive = FALSE)
-#'
-#' demo.obj <- change.clust(demo.obj, clust.reset = TRUE)
-#' cluster.plot(demo.obj,plot.type = "umap",interactive = FALSE)
#' @export
change.clust <- function (x = NULL,
change.clust = 0,
diff --git a/R/F0027.R b/R/F0027.R
index ff900ee..ffaf159 100644
--- a/R/F0027.R
+++ b/R/F0027.R
@@ -14,11 +14,6 @@
#' @param out.name If "interactive" is set to TRUE, the out put name for HTML, default = "plot".
#' @param conds.to.plot Choose the conditions you want to see in the plot, default = NULL (all conditions).
#' @return An object of class iCellR.
-#' @examples
-#' clust.stats.plot(demo.obj,
-#' plot.type = "box.mito",
-#' interactive = FALSE,
-#' out.name = "box.mito.clusters")
#' @export
clust.stats.plot <- function (x = NULL,
plot.type = "box.mito",
diff --git a/R/F0029.R b/R/F0029.R
index 4672c2d..3747386 100644
--- a/R/F0029.R
+++ b/R/F0029.R
@@ -7,9 +7,6 @@
#' @param filt.ambig Filter markers that are seen for more than one cluster, default = TRUE.
#' @param cluster Choose a cluster to find markers for. If 0, it would find markers for all clusters, , default = 0.
#' @return A set of gene names
-#' @examples
-#' marker.genes <- findMarkers(demo.obj,fold.change = 2,padjval = 0.1,uniq = TRUE)
-#' top.markers(marker.genes, topde = 10, min.base.mean = 0.8)
#' @import Matrix
#' @export
top.markers <- function (x = NULL, topde = 10,
diff --git a/R/F0030.R b/R/F0030.R
index d14c398..c462454 100644
--- a/R/F0030.R
+++ b/R/F0030.R
@@ -16,16 +16,6 @@
#' @param cex.col Chhose a size, default = 10.
#' @param cex.row Choose a size, default = 10.
#' @return An object of class iCellR
-#' @examples
-#' marker.genes <- findMarkers(demo.obj,fold.change = 2,padjval = 0.1,uniq = TRUE)
-#'
-#' MyGenes <- top.markers(marker.genes, topde = 10, min.base.mean = 0.8)
-#'
-#' heatmap.gg.plot(demo.obj,
-#' gene = MyGenes,
-#' out.name = "plot",
-#' cluster.by = "clusters",
-#' interactive = FALSE)
#' @import pheatmap
#' @importFrom reshape melt
#' @importFrom htmlwidgets saveWidget
diff --git a/R/F0031.R b/R/F0031.R
index 68cecd2..e246ff9 100644
--- a/R/F0031.R
+++ b/R/F0031.R
@@ -27,21 +27,6 @@
#' @param out.name If "interactive" is set to TRUE, the out put name for HTML, default = "plot".
#' @param write.data Write export the data used for the plot plot, default = TFALSE.
#' @return An object of class iCellR.
-#' @examples
-#' gene.plot(demo.obj, gene = "CD74",interactive = FALSE)
-#'
-#' gene.plot(demo.obj, gene = "CD74",plot.data.type = "umap",interactive = FALSE)
-#'
-#' gene.plot(demo.obj, gene = "CD74",
-#' plot.data.type = "umap",
-#' interactive = FALSE,
-#' plot.type = "barplot")
-#'
-#' gene.plot(demo.obj, gene = "CD74",
-#' plot.data.type = "umap",
-#' interactive = FALSE,
-#' plot.type = "boxplot")
-#'
#' @importFrom ggpubr stat_compare_means
#' @import plyr
#' @importFrom htmlwidgets saveWidget
diff --git a/R/F0032.R b/R/F0032.R
index 706d86f..d3964f8 100644
--- a/R/F0032.R
+++ b/R/F0032.R
@@ -11,16 +11,6 @@
#' @param type Choose from "classic", "jitter", "unrooted", "fan", "cladogram", "radial", default = "classic".
#' @param cex Text size, default = 1.
#' @return An object of class iCellR.
-#' @examples
-#' marker.genes <- findMarkers(demo.obj,fold.change = 2,padjval = 0.1,uniq = TRUE)
-#'
-#' MyGenes <- top.markers(marker.genes, topde = 10, min.base.mean = 0.8)
-#'
-#' pseudotime.tree(demo.obj,
-#' marker.genes = MyGenes,
-#' type = "unrooted",
-#' clust.method = "complete")
-#'
#' @import gridExtra
#' @import ggdendro
#' @import ape
diff --git a/R/F0034.R b/R/F0034.R
index 2050d59..28dd54a 100644
--- a/R/F0034.R
+++ b/R/F0034.R
@@ -10,11 +10,6 @@
#' @param pval.test Choose from "t.test", "wilcox.test", default = "t.test".
#' @param p.adjust.method Correction method. Choose from "holm", "hochberg", "hommel", "bonferroni", "BH", "BY","fdr", "none", default = "hochberg".
#' @return An object of class iCellR
-#' @examples
-#' diff.res <- run.diff.exp(demo.obj, de.by = "clusters", cond.1 = c(1), cond.2 = c(2))
-#'
-#' head(diff.res)
-#'
#' @export
run.diff.exp <- function (x = NULL,
data.type = "main",
diff --git a/R/F0035.R b/R/F0035.R
index 1cbd177..918e632 100644
--- a/R/F0035.R
+++ b/R/F0035.R
@@ -15,22 +15,6 @@
#' @param interactive If set to TRUE an interactive HTML file will be created, default = TRUE.
#' @param out.name If "interactive" is set to TRUE, the output name for HTML, default = "plot".
#' @return Plots
-#' @examples
-#'
-#' diff.res <- run.diff.exp(demo.obj, de.by = "clusters", cond.1 = c(1), cond.2 = c(2))
-#'
-#' volcano.ma.plot(diff.res,
-#' sig.value = "pval",
-#' sig.line = 0.05,
-#' plot.type = "volcano",
-#' interactive = FALSE)
-#'
-#' volcano.ma.plot(diff.res,
-#' sig.value = "pval",
-#' sig.line = 0.05,
-#' plot.type = "ma",
-#' interactive = FALSE)
-#'
#' @importFrom grDevices col2rgb colorRampPalette rgb
#' @importFrom methods new
#' @importFrom stats aggregate as.dendrogram cor cor.test dist hclust p.adjust prcomp quantile sd t.test
diff --git a/R/F0039.R b/R/F0039.R
index f85e3ed..91fa66d 100644
--- a/R/F0039.R
+++ b/R/F0039.R
@@ -4,17 +4,6 @@
#' @param vdj.data A data frame containing vdj information.
#' @param cond.name Conditions.
#' @return An object of class iCellR
-#' @examples
-#' my.vdj <- read.csv(file = system.file('extdata', 'all_contig_annotations.csv',
-#' package = 'iCellR'),
-#' as.is = TRUE)
-#' head(my.vdj)
-#' dim(my.vdj)
-#'
-#' My.VDJ <- prep.vdj(vdj.data = my.vdj, cond.name = "NULL")
-#' head(My.VDJ)
-#' dim(My.VDJ)
-#'
#' @export
prep.vdj <- function (vdj.data = "data.frame", cond.name = "NULL") {
# read VDJ data
diff --git a/R/F0040.R b/R/F0040.R
index b456a7c..b89cde9 100644
--- a/R/F0040.R
+++ b/R/F0040.R
@@ -4,19 +4,6 @@
#' @param x An object of class iCellR.
#' @param vdj.data A data frame containing VDJ information for cells.
#' @return An object of class iCellR
-#' @examples
-#' my.vdj <- read.csv(file = system.file('extdata', 'all_contig_annotations.csv',
-#' package = 'iCellR'),
-#' as.is = TRUE)
-#' head(my.vdj)
-#' dim(my.vdj)
-#'
-#' My.VDJ <- prep.vdj(vdj.data = my.vdj, cond.name = "NULL")
-#' head(My.VDJ)
-#' dim(My.VDJ)
-#'
-#' my.obj <- add.vdj(demo.obj, vdj.data = My.VDJ)
-#'
#' @export
add.vdj <- function (x = NULL, vdj.data = "data.frame") {
if ("iCellR" != class(x)[1]) {
diff --git a/R/F0042.R b/R/F0042.R
index 0fba2e3..894bae9 100644
--- a/R/F0042.R
+++ b/R/F0042.R
@@ -3,19 +3,6 @@
#' This function takes a data frame of VDJ info per cell and dose QC.
#' @param my.vdj A data frame containing VDJ data for cells.
#' @return An object of class iCellR
-#' @examples
-#' my.vdj <- read.csv(file = system.file('extdata', 'all_contig_annotations.csv',
-#' package = 'iCellR'),
-#' as.is = TRUE)
-#' head(my.vdj)
-#' dim(my.vdj)
-#'
-#' My.VDJ <- prep.vdj(vdj.data = my.vdj, cond.name = "NULL")
-#' head(My.VDJ)
-#' dim(My.VDJ)
-#'
-#' vdj.stats(My.VDJ)
-#'
#' @export
vdj.stats <- function (my.vdj = "data.frame") {
# read VDJ data
diff --git a/R/F0043.R b/R/F0043.R
index 07c7473..d48400c 100644
--- a/R/F0043.R
+++ b/R/F0043.R
@@ -24,11 +24,6 @@
#' @param eta numeric; Learning rate (default: 200.0)
#' @param exaggeration_factor numeric; Exaggeration factor used to multiply the P matrix in the first part of the optimization (default: 12.0)
#' @return An object of class iCellR.
-#' @examples
-#' demo.obj <- run.tsne(demo.obj, perplexity = 20)
-#'
-#' head(demo.obj@tsne.data)
-#'
#' @import Rtsne
#' @export
run.tsne <- function (x = NULL,
diff --git a/R/F0044.R b/R/F0044.R
index 2c26641..9ea9274 100644
--- a/R/F0044.R
+++ b/R/F0044.R
@@ -4,6 +4,7 @@ setClass("iCellR", representation (raw.data = "data.frame",
scaled.data = "data.frame",
batch.aligned.data = "data.frame",
metadata = "data.frame",
+ pram = "character",
atac.raw = "data.frame",
atac.main = "data.frame",
atac.imputed = "data.frame",
diff --git a/R/F0052.R b/R/F0052.R
index d572072..ebbc141 100644
--- a/R/F0052.R
+++ b/R/F0052.R
@@ -1,7 +1,26 @@
-#' An object of class iCellR for demo
+#' Make BED Files
#'
-#' A demo object
-#'
-#' @format Subset of the data with 200 genes and 90 cells
-#' @source \url{https://s3-us-west-2.amazonaws.com/10x.files/samples/cell/pbmc3k/pbmc3k_filtered_gene_bc_matrices.tar.gz}
-"demo.obj"
+#' This function takes peak marker files and makes the bed files per cluster.
+#' @param x Peak marker file.
+#' @return Bed files
+#' @export
+make.bed <- function (x = NULL) {
+# get filed
+ results <- subset(x,select=c(gene,clusters,AvExpInOtherClusters))
+ peaks <- as.character(results$gene)
+ peaks <- (gsub("\\.","_",peaks))
+ peaks <- data.frame(do.call('rbind', strsplit(as.character(peaks),'_',fixed=TRUE)))
+ results <- cbind(peaks,results)
+#######
+ My.clusters <- unique(results$clusters)
+#######
+ for(i in My.clusters){
+ dat <- subset(results, results$clusters == i)
+ dat <- subset(dat,select=c(X1,X2,X3,gene,AvExpInOtherClusters))
+ MyName <- paste("peaks_cluster_",i,".bed",sep="")
+ if(dim(dat)[1] > 0) {
+ write.table(dat, MyName, sep="\t", quote = FALSE,row.names =FALSE, col.names = FALSE,)
+ }
+ }
+}
+
diff --git a/R/F0057.R b/R/F0057.R
index 211848d..5c7e585 100644
--- a/R/F0057.R
+++ b/R/F0057.R
@@ -5,17 +5,6 @@
#' @param cov.thr A number which average coverage is divided by to set a threshold for low coverage. For example 10 means it is 10 time less than the average. default = 10.
#' @param assignment.thr A percent above which you decide to set as a good sample assignment/HTO, default = 80.
#' @return An object of class iCellR
-#' @examples
-#' my.hto <- read.table(file = system.file('extdata', 'dense_umis.tsv',
-#' package = 'iCellR'),
-#' as.is = TRUE)
-#' head(my.hto)[1:5]
-#'
-#' htos <- hto.anno(hto.data = my.hto)
-#' head(htos)
-#'
-#' boxplot(htos$percent.match)
-#'
#' @export
hto.anno <- function (hto.data = "data.frame",
cov.thr = 10,
diff --git a/R/F0072.R b/R/F0072.R
deleted file mode 100644
index ebbc141..0000000
--- a/R/F0072.R
+++ /dev/null
@@ -1,26 +0,0 @@
-#' Make BED Files
-#'
-#' This function takes peak marker files and makes the bed files per cluster.
-#' @param x Peak marker file.
-#' @return Bed files
-#' @export
-make.bed <- function (x = NULL) {
-# get filed
- results <- subset(x,select=c(gene,clusters,AvExpInOtherClusters))
- peaks <- as.character(results$gene)
- peaks <- (gsub("\\.","_",peaks))
- peaks <- data.frame(do.call('rbind', strsplit(as.character(peaks),'_',fixed=TRUE)))
- results <- cbind(peaks,results)
-#######
- My.clusters <- unique(results$clusters)
-#######
- for(i in My.clusters){
- dat <- subset(results, results$clusters == i)
- dat <- subset(dat,select=c(X1,X2,X3,gene,AvExpInOtherClusters))
- MyName <- paste("peaks_cluster_",i,".bed",sep="")
- if(dim(dat)[1] > 0) {
- write.table(dat, MyName, sep="\t", quote = FALSE,row.names =FALSE, col.names = FALSE,)
- }
- }
-}
-
diff --git a/R/F00100.R b/R/F0100.R
similarity index 98%
rename from R/F00100.R
rename to R/F0100.R
index 562ea29..678028f 100644
--- a/R/F00100.R
+++ b/R/F0100.R
@@ -156,7 +156,7 @@ utils::globalVariables(c("%>%",
# cell.cycle.R | F0049.R
# clust.ord.R | F0050.R
# data.Sphase.R | F0051.R
-# data.demo.obj.R | F0052.R
+# make.bed.R | F0052.R (Used to be data.demo.obj.R | F0052.R)
# data.g2m.R | F0053.R
# findKNN.R | F0054.R
# gg.cor.R | F0055.R
@@ -176,4 +176,4 @@ utils::globalVariables(c("%>%",
# add.10x.image.R | F0069.R
# spatial.plot.R | F0070.R
# i.score.R | F0071.R
-# make.bed.R | F0072.R
+
diff --git a/data/demo.obj.rda b/data/demo.obj.rda
deleted file mode 100644
index f3d17b0798394e0e658de8df57e313ce3274d836..0000000000000000000000000000000000000000
GIT binary patch
literal 0
HcmV?d00001
literal 257086
zcmagDb8sb2@CW!}+kUZ~jq$>bZEdU<+qP}nw(Si!wry=To6UWH_q)5gf3B`;YWg$X
zJ=NXSHB&XyV!Bp*?4mj}YC1v)>(>ywbOB%ge{@8YFSMuS3<3c5WW^#1J@0eYGK#Jt
zD+>V-AY}l64gi1*{0{>DhyGjt55b}TFTd@5c{u=)fdBxA{6W~@C!4_kd{Ucm6J#-pOX>&~OW$7lT5emA0r
zQGub`&Qi;?Mt!Vn{dRv*z@9MrK^35$m05(_~xpn#qMNb
z`qjR@_3#mqIKJ)H!7z-m`thPKk6+890}Z^U;ePA3T;jmvy1Q;5z_+H**WuBT|JHY8
zvk6Dil~K=J*2Cz)yVjt1yr#1+E96#3!sEepC;MjQ!P3c?-qd%?ShD2NH&Ic$jN<(|
zQQy{}i)MG$o|;k3Cv@l0qt>SpSRe$4sQ&V;dC(z^3wV2`#wdhry!@ZV&`&MCR
z0G@kuSzqhp@DZ{*Z(Ci`>wNP`-U(~8`mYfGORoRJ3<9lS<1{qn!dc3|vLgzG+-v%L
z?b?Xt%hr7r0=v)0TU*{^H$L`Hxo5@KNmq6YJOVj&mj)v`?RM^`D|2%h9)7Og_m%Cp
z6BdkhZr_*PDI!XGG~0K6=GcPUPk(xVdwQIvy0@Bc`mPujsjJMh>L{KGs3jt0(`}MTe0%1_9MkQh
z#nZd%edn{iqWbL9+U3Kz)xK+Y)#sBwH@Df-SE$%_ezmLVVN30%^qg#S#8%cTl)`tC+T^y@N?7Wa;p_DawS8H$bxP)zWqM1U`%OpO((aDB4_IS)*ZYE@XC75xhjDjiMrdc*
zM_|{zXG?!m%|Z5z&BMO^&8Gp}cXR3!AAj4DvzC|Iu`$zcsldmub-QSo^Zk2C`A8xRJ0vx*vwWM~=`rztFU;6?90GR(hS^$W1
z0Kf_j1hHUmXFqp)9;l8^K
zu`*snu0RM0Bw0uE0l1O|J`?rc8Q5(Dp8375p%zG^ix$un%w8_GxB%b^#Xcbn<~9#U
z&6~b#lCw@Y?!`NwAaTRkc-5EHbkn82O_^z?%{OIf^%m_c?0MR{89CyW2UOhb6-G
z#ZqXX830?z8I7AX-hj|&mzV2G^x|r{Jn#0bZXJdnO%}XX(`a2%UsKJ~k$9&K(^8`T
z?A90F(s-*|w*c9+1rh@2L4r^kU6qZk7ld3N7)uKKLJYhS
z1v2XZ&k(z8ZWRhba3H{Chyk2@;4?(g1@)<#DW=+LC6Cc%Hy2!qxb2fw{yI(}kM}m6
zhX+*Q>qJrAh7a!R(;&9c&1#^Ef?SRS64Vzw|0B~FuVT{ua=Wwo$C
zXdU>UB>;f+U(FCvq=Wur7=bj9Y;j!6SRktf6Em)bWlE|4tDpc#D-2&`xCj7<{(*!Lg@_nRkbpoS
zJOBi3T@`FbTNOHOL0{}0cOF+IB9^1}>yx&yA~uS&_MpkQKnroL$qEwA**iz1xX#Qtn
zs$8i_pfaq;^0KcYpK?!mtCEKF+PGe=j+*y#N95F!PZD%#nXCJ^8vKZP2f^GJ;R%)u
z{p_rKIq)yq31Ft1mv4REtbR!a%A(T&_ZC0^xX*2MhD~28_KY23qW6tS=e{$m#yIWgA|zmW
z={PtMJ_|h(G>xDJ_qn#MPzhW_2V!8NM@LehYodz39y$=g|9DoMh*{kt4m^fWQ{FUh
z_%+--$i^%8O4PqUWx?mYTElPJSG;5L(=0x_k2)Wz#q}m(AH9=h%t*H`^rBUb4H4p0
zrkO4T!;RR+lW!(uT|X4$E9wfz*Uib}p`p4n}8)3f8756@x||%<5Bt9*Y0!;gl&@067NJ!*BF4aXM_ZigmF0
z?U$`MLx)N+L1T)N3t^blMpB#g|2#*g0LkL>&JB{|C-p6UFHU42kT(N18sS*=sjF;a
z%`wA3G%?&eHx5C>(b%35h|W>@oIrELBzH>AyfukJ;3KW6WTrd7x@e3;IVL0Mn;B7K
zU@cE-gC!o|G7&o!
zLg90bBx;eB&~ncuGxyt(2IYfJgP_UTA>WQsXQR|pC$glJ*cY
zPga{b9%gJZmK&WIy7rQYr$wKetu7>S2%$;Nu$qtP{5)1;l17RqL2YKsrwb+zFicN0
zt@Itpt8yq|ft=+_yUkXx@Nv8mG>dtW6pJ_{ua5VLay!J-JI@+ugIBMT^;cWOuI{-6
zKQ!dQLcdZK8X4mk2UV9IcszIKS@N&2C$Wja)x$?Cs+OSxL{NpL68rWr^xZJq-+87C
zA(V3xj=t{*I19?(gNTK_3!@`01ghZ?={@6DysNN#p}~Tm(6`9Pm|0%R^fae=c!pJ%&K2rQ#yQhBv-^h(r(Q^1V
zZ&f{R4$}RTIKVQ)>oEwv>e<^N5HbZSx0oz`FF7B>@gw052`3Dbn$E3xv+j@ZR|5i8
z(e=LbYBiL4iDqyP)X-dGi;HFFG?lEhL({t_QxcWSdxk4^OKBJIiVx%Ekhi;en&pdI
zL48n@`nDqf<@hTa-k7`;>5E)9rL{6Sy0Aj~@x@I>fg`8r)~m^~if^=}$>UqbYUG1$
z7p-fJ1j#Xiy~-!P+_y|w`+@vAXXHZ?s)Q<oi<3+5(&Z~cVRir?#Gk1)6TznODl4)`i*&e1uMPhM^J*iO(|y~$GW+`DBR
zBT-%r6Ry6A?JvSUaBi~69z8kk-_h^%WTlXR^dB2K)94nd=H%o$DsATfqFJP6^<*MP
z&@#7cWH#6FM+eyh{iGtqZDV|f%fC$8BcwCSE>)8VFe#?XB~AJi2>5!nw_C)|D>NpV^981H|?2o*~
zYtn;i*3YAi;nM00H&zL80v=VPXk<-bLs2r<)oF6q^n{!{Ni%BU->b&3V)?(KmlMqi
zO^Ge^3o$bavR)nPMDahMY+1L2+Okd_O7FXN2ZUj5jKC0U%Bqc*L>bock;1xM}?WJ
zB3N8{s_R}E$=0k;6<6SCSo@GR9%huS323M>TN9V--XEE-@$T?jxhEau3FtCL50)v8
zHK=GAIB9y^$Y%43YFrM&mYdYl*Z^fM0EJT0(Nv-3ur5Y4swR}wI9`>E0oeQdK2*p^
zzYnEfKta5sq3F=E(w3$eNifjpq$k*!1Q3yEsZps@r$OJ^qO0}1zQW!H)J5#>`Yjt)
z&L__PCZPS;{p-BdV+|I9^M96umH=f0-f+PYl3Y?6Q+cSXo4KJ;XTBH-B5~v~`fB!TnWSAyeq9Bmo
znr75pY%np%AGR03S!dhF6#YzJ{Kw^7ID+Z;^sJ#V
zB0Ya&_9gg(9XiQSlDRpgQ<*dxx`f}9Vm}U5mi>~krrLzWdsx0WBejtXJ3j0*9f#i8
zNKLYBjN=F6rA{d}XEu&|+Q{2dkZ2r)aa1HbZ4*YQ5N{KMR0;Ksc2<4jCSv4%vo5bO
z%U!Yh9>#o1F;Oz-a@ltvk4VVaPejB`@eR>n6;~Q{>IHti_R044i#_Zu&kL4m^kLqk
z{Ebqw3oh1b$5`Te`qOWv-^^S&9f(Xp6s)w=#oV;~Ll!H0nJT*juQt)yiN!NawMssWR!Dw4Ite-bWVihK{4<7PTm(Kq;HL)~M-)BKAi|^RXFBFEvqD
zNdU{($u7_qM_Dk{zRI$Rv^PLyZT548gmMz{Mtb6v9^yvy+a+a)lB4QTmUtbn=9ZqJ
zt{Q<28d-ulSAYave%pDjf>~Y#pMm5bvkatM_l$gl#|9L!4$2EBtZO=nZRhC#F@uVenES
zFozEu5F$pi1p)^q)`<1S)uywMj#I(#RVjFc^@Wg^|e@1(GhWgrqQ0CFn8?
zEe)>CuLcxVcRdYpndM`z;zYecav7n1
zb$n(@ilIU_L4&5})2O+s+;aT|(UXDf2(kEV;)HUDmpy0^DEUBWNV1#2Kwo3HLM*{T
zX;W*qY;0)E(BN7NoFM+I6)@10K27#LVI8T?zkz=+$W6QWGUHK=4TBx{jJ|0LPXQD+
zP52#&qb!bDnKT=0O=I;9e{Kva*v649l0iNW22!@XRb>1-@lYw9>9?#vfy6X@cg)Uz
z$1`uiVkLN?ij!#1p@@U7@*`2OqoR?$aL#PG6Gh%d;8p5c27IYaXILSY&hEx6&mMwQ
zWq63gw?bf9vYL1HIz>bhg?6bJoY7Q%IR$wm|Kc1BI}OM=Slj76yyEYe&s1>9m83?b
zjPh+%b=v%Q#vdH=_zjBNf0(V|`@6>xL!`(#T4%L&R2R9Z)yXNv@L!4wHACSd1DVh#
zE+ZjEQ|V}2wX1c5prD;1Gw~6M3OA*I(XtY2L<8B{orUO_oOd&MTCI_~ffOKCa<;RY
zN9mGR3x>qndj?2`NYwD`wY(?(Zy4v}?QD~}YOJRi6joFw@HD37OHBaFKjHS7L30(s
zkj;|}^)XpdvyKvSk_JVpowfVyg}KVLZ=PM4lH~de^RKb6CTw(=6bKkx9aG}NoB7a$26&L;T7v92g7G_vp2mOJ1Tb>d
zHGm(@-z+Qcnt*-TCy}v;p!+_8NS(Mt@?n}v8atT5J#uj~D0J3ju!XZ$t>Ep?s1pVR
zv^&JTo0sf4Sh5<>zWJZOZzTu3nMg}L+h}8hcT)e(1;;5U?a<6%P`Rpl^=-1E!t7IW
zTU{xt19N=X;(pMfN6ZCWN~cNmM0+9=lSZTj)Q$~5_y{$9S9Wp@1;_op_@Z1~BnR>Z
z4rDQX^1^1{{=2DvwGmbT(Ux-MSXq@2<$swEPV}jXL6z54vVi)Dgi2JV?!)8{mL(X{L`7!
z2|95iKo?I8{FKrzI+#2+RO2b=@Y1QYbE(Y762+$lO=81i2Nr@L0gx2j+QIb-}odYP`xp~K%d50I~K-BQ-bsn
z#ygpz&wx08p@Cwd6mp{Y8z)7Uyz%-+Q>xz?~MJ5tG?ord
zx|k>qK9{Nqgd@+4%MAyEPIjL-jC??j-3%&ABEmwpumAwbO(X&(;0>kW%+L{(2Zs(V
zV@GtYoKb#c&OGLt4poK4Bp=`Pl%c+?C>Ri_mzqQst>_zjHD>QdxzxL+hK2gc)UP%k
zp9?B)qv${|N7{rRBEWNF)#%XjOKT#dUq@ow3dqc1Qtuvj8jgw^nMlxJ^5Ij5Ulam?
z1n3Bo=(vZ!?iMVNF(F|^2xRspAjzFH$|oQ4;K`&172vOZw||!e>c~Y^3-wzU`qcTt
z96LDCw(=D*ENard!XZ{BX3_l>xb-PM>@=33i6&o9K&Mh;v!t=zl}hS$LZ?lF|0m;C
z$a_k|L_t9oX(?KV983I`$&hvp-3fv@ot725A`Xvz=At90#^;m74)G;{Yr!Nl?n4~Y
zC%YKZ2)d)MM?x&ogVKhmP1Awjf$%5;9!CWMa$6vi7Xw3H=rQ07ZID&zA=c9=>7$&Z
zkXDp24U&R&4|xQa&_xp9bd{Bz6Y)nDsFx(@ZE1aRMlrlh;yEDOfrAc=S-Io3Hrwb)
zo9b=jG?<}uQo$-tQ0$F0aEBGhkGVzgm|qf2G5K>okeOW@)Ya9u?*)}zzu8I(Y4)-p
ziJ)zDn#_%V%A?@GiIWq?FcNIIvXm6YnF*TVM2}60um8d&&rsqD)NW!e`TItvE@hjR
zT#fi0)*u0n9>oMT*2w*SaFnyb<1H<6Nq6A%u)R~?R4Fl1B)KT51XZHU{xVOCnG*1b
z0Uz6D1O9ShdF^{d(ioxnobPr0`2!97`~B06lE7x9$tODC9rJenWQ{bgNVuwm%B>Hk
zHSjv%nS=6^GSn~r^1@b(65?!sU~b^WMv38Rj7=$yj7{YG1%MNX%cywm(9Qs$AdVn>
zGuJ>q-oK8dTOj$%fJ1mTvR$zYV)IKW{8;)U{IMsxzO~|xg7JoN+I~4D3Ttx)X*loy
zDueWT2D8nQ8gZNJ-|U3_7BJMc8fq@B;}mUZJ1MRLjl^UE-HeSF
z0tbab>B+>7pOq|^w;qNrXs>L+7#S8LFuTd>_&wcwx0n{HR_2h3OS2Q^cV}AEDRp5W
zr)#rG5N0kq_P~j1kdOaPgKY_R7ETB(IR-~}cregBrn?M}oNqpagCv|Bd3MP{@jC@-
zgE!r(PdEMIzA7Fdb2QMp>ru)P9;=u_t{TJVoog6x+6dJP9b1}TZNKKgS6H6#92U|9
z$m+ION=)(X8MIDkMtKGlP3%8}#)V{PCSpj$r<;ldiNj$w(v^9+bsqrmI*H+-9XbmI
z@X6YoyrtGCFrgbce|2Qa+`0CS23d(IJ5e;oj4x3_v0-6CVaSp*1orZ82Uvy0Fs5RI
zMyaG^pp>zra3UDLHrWF22#}^u5}b)MD@vb7`Sg*r^~-E(IO
zKOK~3_fvtc$3FM&Z*~S1(E2c-ikUv7t&NS5S;(>^RB8*|?EC2E)cLH;eUR3hD$1tX
zz2Ogjwh5?Xr$k-NpV>hBL?ZP{B*nvaKKe~yojA1|dZ_Kp<>4N|VT(`rK?JRIb_Qqv
zRGU3x{^XTn|9ps;F!$8Qo=0YHI}8OfF@#w^ybPH3%
zQfzbaRz&W(c?J}e-ytbpEF=?M{7qP(^AM~!&+{^bFeG^3Q*q=*b?G;DIArFL*3yfTN7
zr7b_qJOlWyPf#yX!{a!;X)MBJ5b5a1;TclaxU4XNiR@4bzKd(Ni53pcB&j@CYUYef
zmDKH$48p(ha(BxTl@XYQ$0T$7wb
z=$dg8eD7ik|4sIu7XDJSvoUbolpZqCpTTRx+6Me?QD%
zc9U%O$s?AKP84{AY187uh#qXo6ogiNiox8B@NCH!nh@B-gZ?I?Dmkr1-EGqMnXT&M
zFT_|)b#?W(O2eS)w~^Plt>RxTtVT`ed;0L6!UK{4V67>U&luGwSZ2mcM0=)qGk44(
zIGc+{$h#~dBrH_xtd@&+H(5W2P
z`Ei0b%4Gemh`m!y{-j1-bbGp64{Q_Kt9s}mMfqTYP{VZQPGpL9(&9S9y)WF6f|hiE
zKOEC&XY!9@6CHsZmYpgP!3X~ihyC=&WEnOs?Aj>wCVH?q9!Q=j*Y4XQyi(|9GX#HZ
zHPy?pGh|E_vs`S5Dt!Z#0KuDY7!;X!pU0353-kziukPveVT73&!jmJ&Gg**OzTLSQ
zS!R$%_|pnczn(dGB0@=E!{^j{oc_mCGm@LOUd
zTA4zpmNdrVq%b1Hgs2qCFO8%@R&3qI#}qHS_oa)RfD5bq3*SX&tc_zOlIw$edFVNffHr&m@a$@Ej(G+hNPxSaq4VL&cZECKgSw
zBTL35kYo)#mt|Ok%$B*UOZ1>dWV#Rz4Mh&GVv~S^p*&NNvmS4sg|4z&c+ic{dfCJe
zGYu^*GpyQ)mGyC~?@_#~xe|x#+PhbQ?U>%i{9bz(UJ(TjQ}8rs-(UlFhe6Yw8+D@n
zNnSWX+UEgMu=Js(G?_OC!P{XmLl$3klDQRVf;
z@lg1qaS&fG@9GcMipSahUMYI(+$JJjb13F0@#8js00dhDIf%s6G9o(x6_q;}5swhQ
zV2okvKS>??2sK(-O1UHo5TK%XhM$xX>l@>uL(8qMhAlm@NMg*hprk<{p6oGvI?p|!
z70@^(XC^R7{;p4LjhMQSgW0&aBetYX6p?i4Vuo&OalnCxf2h^4A#%uILeER6H|Z%9
zoPNcz{5_5saeDs6@r2uptzgY&hu~0*^)F=IyJ}*gKxua;5P$5QXcbzb8s1O@xy?oD
zvxIFt?<)H&PJWm5CTJ8JA%cl7uU$E3n4c*h2KOmIAc?V7z(FqBZRmBnSfWTP^)}twW+A7|x
z4>ap3tN+Fx?m^JI8c(mBmEy9ksB6rU+5&$wSb9w@NlGv}c1A$577`IM<(w2|V5Nn)
zVD_eaSQ!${gh#Gd^#b`;I?0r{Ce7wX*rNNQ`u&sf
z#!G_Ag(}1p?k~raC+}gN>o&oKTl>Tf#e)JNu@S}&AbtDR4U)QyQp*VoOTd+H#P@ab
z070$!O>g`k-kH=W%bB;CYRvM!)JD2C4&GKPSB8W=`~4}lK?m(W8icdt1U}8s&J>HXrR$|rM(whHxHDil?zA
zLP8nJ*y%fu%uP!+ZnkF9N);&G3x%J99+a6pGRijn6l22SyM-3~Zl2RowXZ>h>{v#_
zaUV`OveUJ_!GE4HND_}=#GEjld9To
z*ZZT08Td(Al6T@uU;k%;BObx$+&e)zj>GM{ny?!s$K5Xj;?65Z>!bM9e?hq^FAT
z(Ve@}tW+L{Swjh5Od~x4KSS5Y^oS({7gdfRU2jW?VB4S@Do|6+gYLW}!==8imP9ci
zRrfz6?il#C>e_W#Dq-FVuQMc2i`-x@*soaF`=Ydg37LNBwM$p;`vlc9s=Z!yu)_b4
zUya6jB81BRYcse42Tv!M&Zkjp5m<2<&_NFBJU`!qt7VhuFdzBOP~Obo;_Pp`bEq8P
z!kP~*g_2^_DFaH>AIMu!4g9|>Fo;Jhxvt0GQlQ&EUi$sIBp8}($9FnxH*XSG=nuC&w@_vGhR@}eg^mi?+$gmQv>y6tRBor*ck+mz)+X})DR`H$>^b~s
zKqcbE6Q&T*D}D39y>KtId@#s*(sI5nKYI#S2Q=R1
z)0%56HnswJUJ*0;nr@>cXgS&6MpnNR-eyxQ!#vymR2Wa
z{&CztQ!F-5_2Ly9`Jq&KKSy{|{%wYEs{#DUvPOz|T+I2_x#k864Hm*<$vk2fuix#fz5nO;t%7g;25H}v-E5~*L22ic&Ee3eEA^V9t!t1Xxjnvz^-RzF
zen9=8HXIV~C(2D#j+41H|jvE5hL-cqDj6#;CQ{p(A?FANMHF+2}uXuu$zS9q@u3s6VyiJJjI5b2D8F`FxhmSl}KZbu}E=hTkX+^P-^|_5}H&
z^-?bhLB7uJULIX}Ap#
zlnbbg+8S}-vFk%(ZlCkcG15CIvzx$opV`;O)y6df=+XW2a^Rt0q!Fq
z2xZh$bnLzWw~M=UWm4hs{V5`kN7g$cKKu=;Pn+7Qrp;(_7r
zS$#2`cc_`Fh~}9(v(jT#e0LmaX_;95kn4y%tKYx$$YpAT97NP6JJ>KHgjXS3)m^V3
zMUB`@h*(fVjqCBIB{tN2BD*p-QnDJ^q7DBxv$hdH=Lv~x8evxK0tL*IgNy)48cSF?
zNVh`L{cfl^Z6%r4CWVnnoQ**XP73)zBdy5|RaBlZpqG*yGcD-c#_3IAMXt`*=zJJ%
zvEznio%&OOB-pKq-Hn}}gEu6^bLRPRApM8%{6s<7U$ArC;o<6ezJt_YaW4C5XCR$1
zW}Y_x&N^r2XXxXM+85)-bnohYMz^5&x$T9M+Y}e`&*_qe)Lx63pLZVb?(+&S$LVi>
z4^T4hogoJ0k@Ijg6vDj|gCxQRwR54-bpmR>$6%!)D{F59jFf`NVV6UG?NhD`e-ODX
z-$by-2JE2sw_vGRV$3XEdAB6gG`7OcVG(nYry(hj(2}JLcC=`!6Q|2XEk0_s{r7DL
zsZZA4A^b_>)^?@n=R_=+ZJgoTa&b1o_h`Gy9h)xB=I_Hk&2KX!#HA_8yZgwBNu{OX
z5BTiYud45I;GHvO{=SnG&IgNPVdG(~>sxt~QE=<+X*sE7;5JZyh_>2seK1v><01dF
zwnYdxMASLwoBx4URzY3Y05fm6r>DptmU|>Hy~2{Q;V>ar(7tW`c;D`eJdc``GgyYv
z!M*-`vT;e$)_qd6;r`6v$d_J_AKhMiIC8`Rz1mHaBkk4ueAk&WrT0(YKm7&pU;8JGbN4SS_87{wU5{iqcN~F|YM97~)Ppu@z~EoNy9nd`%!xl;rbmOr2%JCd
zHrCqRc<}kG;p-pqXCZ3JkR6l$*5f*LF
zfvVY$&9%%(&W8Ni_ChbA@nk{!Q>cCTQr$T3jWGsuKmQf)veYdp1M#r-+MRKF%QCR8
z1yWv0<$T7XYFF%w$)OYze27gf@YpfIF#VUcfFN(eMQC-u?|wQ%@>dzx{MhCyF82mU
z&~5a^>(QeZI9X1l`qCAANTArc*Vg{-uzg)ava7W1^2@fR7n*F4vhl@5lm38sw8O^9
z3UWBGmRMg(OZtUqq9IxT#Zsu(sq5hY+*SEZm0|I7ttWRUaj2%;eU!v+5y{+L?M@1{
z@BTC;);#RMd(AH$zM~f_;CcmXTNt}^XWgC@_w0m;DId4`
zL84ss#pQ4wja>p8{jd+E`eO1A_oWQacDz-8`;`bU?>9sBwFqVh?GF3CCQ>;co&nGB
zYjBc*0P+X7^-zFE+~1PLzGpYioBUb5;)FeqD)+_C)?->~^n&X|mI1Q(x1mw)I9BnE
z>4)}%gR>5O%gDaf$=rWM<1l654(B%~hMzBJ)5$hcf`6lsTHM|C6BU)0i0@C%nLPN`
zr7BPEH$3c#Dl&~5i2T14cjvj$j&&S(2MPk}0T9iDfJOl1jqIoI=|86(xR`>yHN>8j
zNfn=>K~2w<`LrD^!SWB)kiqwbS5csrjj%R!$w0D@wr1uDnr?i04D8_g1FZZ-UM23q<@_&BqX1t*W3aqke+*DrJGZ~
zgfc7VuViR%?b3E;sJi|+A5qiyW~|%Ae%i73@5L!E2{q<^wlQtXCYO~tz0mo{NE;Ub
z3Vn%$~<_z>#Zu=7yj{5&PtXgwAYDFY(4Rbxb$~<2$
zZ@=hOt@_X5+o+yg=(%fYi)(2b5VM{f-;c5?67yXQd1x}F0TbV6Pin=}?t27#IHS*q
z<9EO<1!Q_=uDM0l@tu3MoXw<-ZF|8>8JgiFUOTluON%&%Hra8X2&9(HP#5e@Ud@l5
zettcYABMV!6FQG~EJi5uj%^itFXhy+-%F$
zs7UHO^}kWLk|32SH!O%=+UPy;Qn{7#&oZP>f=jxRKzW{~9CwZ#13M`;kB^S18*8W3
zA09tgK}QYV8xh-%fabw-|DLPbgPlJN9&iFK&=ENiXw(d(q%tv7oyJPcVdKgc2AcDv
z@^HC3EE%RQdB0SKFY(~Eejd!MXPM`EDz$(3@Ou~+rq7i!Y4z0OiL*Y}yiHMEgU9aHsUXK9O@(nNMh$ag)
zPcX;fv*Jy@NXWwT!}m!)10(>OhdXCy
z0~DlL?)}7jus@#-I_B4(Y<5_L-%3AUq8ISvb$d=cOkb#yhv$Sp&xQ%`ENk;o-VS$o
zlpM6Ub$!zDYJ4kG3X*q|UVcA$*+!5mZULU_^8QWT-}7}k|9*25&;off<^LQI&UEtv
z-UzN7ftwz*g@S)Iespu&&TZG9G45SHiv#OhY5G{gSQRjrAmwNj}
zBy&G(kaN%6r(6?6G%|FT&eyQwb5{*U8P4BcAFj_|Z3_(ZXDw^1Zk|3?@gzxhcB)~H
zH$3{x1pLsMeAeF&t~LS+yv%<-*rwXqJJi3uIx0K&9Bgcj3pm|8{M8e%%=6
z(1gA0r{xyIXJ)^bA?NhUDkEhN*OI%vn_jf;DXzYd;QA5>KGQ==oJ}_pG%SXWPAX^{
zr9{6FMMOdB@76}4uY++LS-1srLf#5L%)MF^eV2#LrNt4Rohp=rdo~`OP<*Z*N9GD+
zqN@Ijh=T`Fw$J>#cWSO7TM6O%#9s^xt0;Z~H=IpE7ZZ(oZ!x3b%Xqxvidvj)yrs1p
zo|)(3TD}pu?!*jSRO7f><^G@{xF~O(`^7+J7^?zO<){@RK{h`w9_U3pUi$?4THcAl
z-?>n8lB)kCYd|L4u;+^f6e)iD=;f(zq$+f*5&KHe(@y+s!0B~tN2PdNK`W6lbRds6XCKeW_8ipt09
zR0*vkFjVlz3MUb?yZ`Rb0=&Pcb?G|B{v*|Ykn?oYj(rWDo4be@YTv{!goa09k3xX$
z9HI);9Q(}IW)VIBVq&;YKs%FZIRDqT4YYZ0b@*$ElxpBw=+n3DJXq6D)iUP?Yx^q?VHRj;dm!mDxDTp@Fz!OVFuLCg^|=`~@IXEIt69hKFXx
zYOr}aAO%;G<+moTFbylT_6p94mxo-rPV)KMwGPUIaoBPU8z*2z9@J7~A_)_TOROz6
z8a>=;5ZOFraoEckW(*28nHY=4MX#k}Je(L48$pPf)?UFy9H|u>)<%LbLxM4?#R+cI
zL>$s){2z!!ZpIXr*o-SSsq_B<+5R8f|JSzv7lb)LLS}+6{-4PZ7tnw>+=+9b&lncS
z80N@;gTy@%mO_Hcq%Agrh=k3gT#Gn7ia1*3UR>(y0|D7T6Q0tgKGQr^lqW)p
znSzorIXP1(Fadf2Ef5Hd+0axO+e6@58H&tM_-}p=%BAoE$3x_WMo{R9qcq7NC>fei
z1e#=ga;k+C2=Yh;ctX$;e&Rw6NRJ{bf}
z)$UQ(8cO98A}xk^J_&@dMV5irE*T2Ogeoh!j4Xf-
zM9nING&q7l$h4q&GnW>Fnm=aU*9N704X8*<{|I^>4I+@dWm}C3ah8@kbUzk@qH9GD
z5ubR*41ydp9%R1+G0y|;WYIIImN25;2mU06T$sQF{1o+0OuC%kb4vPck>p2JI|OUH
z#9k}L7)lxhV(Kr3Df5w{*b4Ath5Vsq2aT__PZ^d(Za7zrQN^W=1>GhP36U`w>L9dp
zV|7rEZg}wu9)_Gm;)=%vyFxfjRFJ%d8AJeu_>k8FoU}KAK#F
z8d-zwG9fH+4gao7`r85;jlb3H22ppso&K-+O?4Dj7ns(Y-Vi68bJG}5FlZLlf{}lX=w;TYJ#*DH>^6dMkummOF_O12n
zRUYk$7;#Nz^U)T#+ja|b~>3vclcj_VppUQap;+6~`uY
zGkt8=msF!{Ar%L$rXPYfW3aEXXUd{0m922*j{8%NAm9~SgKRTZ2Is;Cmp`O$*NTL>kCzzx$Xsq&?(4$
ztis~t$oi@9$3p}j7l+Z<#={?PmgYy^{^WEhoy9WA10xAM!^sORDpazxTmBfw_%R|L
zTye!Om|qP}{>HMK2mKkf?Rt@>;m&?zy_nNx_p!@|Xj6BR;Ta!RVPx3s$V=XpD$W?`
zt|u{iX~mAW)ETG0RM&L6Ro<27D(Ng{AcGK=l+)JeL`$1H=@i_sMPLx(l0}d&2Th1W
z`N#Kg8B2txxwZ!oQ~D72s=hN>oS0oVcIH#yOpVnY*nNjI>}UP
z|LOIlbaL%T=&uvOG;kdk45icPq7y2Yc=Y2X&q6WHiFL@sD^SItTfesi~s)$LHH@0uJ)7s~%1~
zZ>!Xy(`l^Um-8v_ha&hdsgiL#jUsSgJpYlaHD{z>6vyRmi6p_dYNQyqS%BAb#{POW
z3o2=q6oZ@S{B**W*l@u4CwljnV<4ZxjES(5@Fj|7<2sLU!<{IdqJ$62(CS!6-)
z!>_eOi(`Z{Um%v_2yWDj(=xzi_RAh)`CLK2P0
zzn~0jCL4R3E58^0TF05yOC$3Vt%63oT+CJza33PbLsP9VHCd$b>5l%?@a1LbOJLvk
z28DEjC0UuAOdjOXB|kmOyF0VvP)AX$_{+X|~Y|RSmj<3gxIoYQ>;a2*I|~YMgB`
z5s}75CeYY6*gl=Xu@u5I5wVS=HzZBisKzG7!BEs-*ftc!RS+?QY+_?7M$k5!K}=T?
z2&SV++7W?HNHoA>ZAyyFRvKlkQIy8u$Tv|DrS3|m(v1y;q(@BeGMi33wbN|5#*=C?
z(kfyb6px;Y}85
zG=@nWk=~Is+LZ#TAx78@q}bREwiM%RrblxlLuhu8VA$FQ$!L*GZ84_E#Z4yBQVF#h
zLtqV&qM9k7Y$Hjq2;)PQQJ{_kvr1GT#1?UinP*Hzm_{}^5ZaO=)dd=}MKoGPW{_;!L$vHh-`O61Z}o95u!3Q+8atSHrs8DgJ77!x@2P-Z6?OZ
z)MH|%Mw=Bdjg1Y8%_zpk+emCFf^2RyMx#S&Fe3=#0yGhlL^hLlWEx4gXeiqoNwlog
zZKgEaVB2CNL}O|&8&GY7X|arKY`Q27iL?&TL28YrG@DVNjkGp4A&desxSDKET}3b(
z8e$4)QHsq4Fs74W+i4(ev9#EmNR1J;Ich10rV*su8WUn|HW8@7w;(p7Y#T{7lWmA>
zN-!gB4WSLTAQ;dlG$ztbp|qIN8*N4rq}Vo_7}(PoAYU1&=VbwsyxKt=l|U{};*{52
zXqu%td@X9$sLvQEZAzZ__LQdP9WgblI?AQtH$gO&li4wp!bOKLeq=_cOuC2L?r&X|
zJwU-92F8->rf=B2?<&nB+mx^8(XIQ7u#{A+CsYIy#*8+>|V&BCWrlKd1=kitH|h8QlkqN2V1XCZ(1gJWHVi!GOsa~n2k%i-iPj)8YO1Y&
z)>W}=Gj&T^xRUtpDzntFv)qT(pgy$JR6qM1X-9>ccB2;3-zl
zG!>75@QZKw0-Jx%1afyd>X!I3=rl6AJUv6<-#wtbg3$HE&fIGvu_}>saqu;
zrW&@p;ITn3hRIfpIotIoZOVD7|66*~9-Uh_F+}Pe;{sP{#j=AV2D#aT7XfTlhD7S0
zs#VCc?i7KR=z_}lz|Hg6GDSyvIq~m{H=nHVtr+2{&Ts6M!EJw=_CAZs!IWCBVW-c4
zH`KU>?7Wi2f3H7O+0M(A>i*vs6NchDe;DK6P(SBW2K2duEXi~utC?#9Zqd~_L#HFP
zixpV9r&ffjE~#m19TiQ;+$u)O8>(gw51GTmYd`kx29OIEoMLvL$1TG1F*(1hZ%h1L
zw-%qwd2U1gw2R%1?PD2;7}=w+jW*k04d=}$FixCJp|B)_*5ZR!Wl&m>Mh?;pY#Z;`
zh}%uJtq2XK+d$hI6liRXuxuObWJ4I(*-|02+ij-N*9O)wQko128%?&8X|y)jqZA#?
z7~^eVG!cbhRXfPF7I$5BUy!Rt!CzraFSgu+Xl*tPH5gMxG@DFfYy{dHE@h~Vp|u)}
z9i(V%O_dT+q8)GC%ux%$Gqik&=IKjU<@CwAeN#!Lhc)7)~^k
zVokBO3~hrNLt>@@7!9OGf(;3@#*&vnBN0s|{C1RJ+ijr@h#PDZVi?3WBsMyTgBwkT
z0tDES29ARbi}dkv2txd0!N=4KY}lwE3(G1>f9mD;q$8E&CSS%@GKiuycdN&(%ZhUp
zDmNP_@n91;z!2Q?Zxzj(I(rR&O;55~*@)e$s}=GLQztB>IGxm3$yRj2`d=3p1QtRf
zt^Q-<5MhQOgOO(IjbL`(!$bA8*Qt%D%@gbmMCkoOekQQwHQ{67;{46JAe@?LM{qriv_TXbV^^Jrl{L>Zx_2&9lHTY(G*JssXVs^aoAlY(kgD|j#
zJct~#3Sm=^rS1Eebxkbgf)(nj8xRU45#|yZzFX0WKoeD>waBzlAbjp>D#ZIci#RXa
zDpn;y7i(9+mS|*D1z&tM#)KfchVP8@IUyf|%2kzdW!XBJC19)cg+U0HH!fcPNh#|x
zW#J4gEj=5PL@VW`7EFzUP&J~Tbs>_=`dTaoyC-8cTN58g*nc9<`L~feqeIn1mK`GQ
zTr!IFT{D+j4+!QDv~}uYlYODuf{EiZhR(
zbAJX*W#IEEyM%miLEjq1@7#)TusYNC8{=
zlCk1sHV({DV)025<`1yZ3}a|GjnK5sU5d7GA_AW!c%5KQr(}B+AAVf1u^DDr*niw!
zm~8lJ-8J(kwksH3hHOLc>Y~AY;pV01HpbY&ox@{rY>iZs$6J(4LpFU%7~j9V^#1dL
z&yV$i{EsNX6B1!qYKTGugdidrCKyC!)$e)rIW%`|marVK4q^koppyoB1s3f80<{sa
z_?SVqHk4#yZ8n>4GlUbUf-;#JY$Ste6w_cxKp?>18BHLjHj{CyFt3gkQ4DAX#x-dZ
zauP^3A+#Rsjgm&&3(T`XQ$uKWiAGawx}pKL(ApaTv^IiW2$&6|sg%@d
zQwGM+HiS76jf@*W>L@X>4UM5jMkdfk1Y~GLXk!&ffs6*yMKBdGG{zTNC372NO{Cb^
z*nzP&lWH)GZLzRJVvLQkw1%5>(4#;$Hj->bG)9Ke8bKNc&^DV%wAc-#qg2vGI$;INXkz@yYFkF5!d!RO}=ndsi%HgSsweKFabT%5DL2lE~-~AG>mp&Mx
zsO9ySbr}n=&k72A4(%OYu@ca5xHJl|J8r!LfA=Lc9xfwdxmj`)JjEs?N(@
znX$0+?L~oq@))y*`$(1=wUI>m|y*o0Z
z-Fl;(4~X00n8@`@gLQu{ql>q{!dcJ9`cU46f``PM+sV^9U06q#M(YN^o;{hP@n$!3
zJk%4iPr-%}G%h-b_{I1vp0K^2mKqz&*v-oDrni#OP>kL&ym((btwpL;F&kXdj
z);K|Rlks4G9xQ;HtZNfUn(S(jV15N{uEzg=6M^)d5{eO&<$xlALWwJ?c3{(pSv4;I
za$L)fXP__4Y46JC6uiWuFxb{?!kN@yfGcVfs|Z705YPyrFhFPpQPUTpNHLG_hxxM4
z#&lNVC8qKDBcP7@gV->c1rh7vq{QUx2H|rBVC}=?`(r#fzL%;mb-I?nhMJQaJvla(
zN{UGxI&|3;RXUBhq>vpiU3TjrKru-^#yo7bCEu3r$(!-XGPlLoiuxqJDz#uvGkE->;@Oj0
zi-y|!>$ix=r56`$E8pwuRk4Vip#xYBJfX+LJFhS-qUnApTN
z*ftTMjg2PKIi_k(tTH4vHrga0Mw4hGEXcvOz^0>ZCmcv@jHV(q4Tz?j2Ec5CLW~;{
zZLw1b#=vbhklGsrn?qpQO{Cfc*oMKh$4w~0NZ5wsjBIQ}VA^br5uvo&NwkF+L5Wk1
z2vouh0ByF>PBa5dt5HfZMhUcvX*Q7?X(r1sDTOqFw4)IilVT9sLTwGUl1;XOv;$x^
z3A72ZHyRXTTL)PZR=OhQAq|a!8$t%tVr?WO2sFL?({&QaGX^8qVvQlBG##@-DWpW@z4o9Eg;KHe5ULQ5?kbs4brRn6u8H31Xc->VosiIm_YC^>^)k?KNgR_HVaD}AAYj^Y9H=EF#F~*#ZJYP0L
zgE~Wo_s{TtdjaVlNyTQ`I`95Y3pC(&&ud%jJ52hFb4CWUzU6jHt!Wv{L6oo#U7;LH
zNbMRzTb{Y2QWyAQ7iiW@!7aILz9y1X(n!m^I4hV4;dL*9foTjA|GFIs0yK#b8%$Q{
z&b`905YltIQCLB*)?^s8A!0}skb#gRbSaFz*jquWn*9hkZDtFy!9Ge5c1%}uQ6vd>
zy@8eaejXA*y@X^SLI#rR0zN(J!4zrN5Zzr;D-JisN2y{=!^qtV1mFC5CvmAXdj*Ksn#|b%A?MlCGP^~ist$ON$NkY5F#y#R?C7wz08%Ga
z(?8Q)-`_N7^G;HVpJ{0DQoGb^nRY4bn>T+vsmF{X97H23m@?d3Rc`pG2t&bDY~OV|
zrNm?45^@rXI&&FuosDQy1!5Kd2|a0$-(l6e>*6c}vwbt81kTuBu%f-z*(Gg8bwRVx
zqr4_zh|!bzZ2=_?<=>{Jtg?e6vOIxOHnjP5FFQKMhLLqWIB$o4*F*F7)Lq0PVp5g5#-y`eCX4kTDb3r
zZm)+c*^kX%9bjjA5L7vKZJS>Ef8!}rZ3{vy4hR6Y1wB@D^xXi?pZZV~+)29k5$q(3
zHNondW#NP&K~~ZxKKN9sMV1$~c1kZ^kaCZUgN`^+VNfT(?rs-@4SC
z++R)VnfuI4(2Rr!C6)w;v!FsW8|qr3QXoi_q)5=%*fz(zF%V$JklIayZGjj`wB4FWj0DbcS$G@O64!^h+-&i(g
zX1joG)DzC*n!(VLC`Bl6&=Nah;@E~@J{O%z@hhV!tUN^4GloJvpXP#M;wz7HLy)Wgx~28bk$|NPiAm$Lt-xg-Z=H5Jzi1~`V(IjDN`
z6!+%us{g;^fDjB|yO>ae(a^Yh7u)lspe)q^wipyG5sJt&uDXPk2@R}Di1H2t^$IL)
zM$48KU?bE+4gC)%N*RcR0t4z$NMj}}ARIiBidz;C4cvXI{9?V{VR&RAONsz~-(C|W
z)-^gn7y2{L>St5R{omsD(e3USMk8=bAk=>Bn{&1~ECzkULmK`@5sA*4Kp3*O3a(5p
zGup%Bl|LUc^nKV;w+O7~s55QS+9a_waj+%g)@Owkh#a_|(Ws5N
z^5aCZ@PUD2*l01VZ%g5ifa4q+0Ld0+N$L{r?Ab7xpKohDwfC)uaghM35LXksC?V&;
zS{k(S>jim->aulZSNLh1O*ENlEhrFWMfXR6+g#b+zlXOyS0+ZvI|j&l{I4TWIKG)?
z%M32}9J;=maPIzjn(U;}Kq&i<+yCY5K5A4_YJq#9cA>AT?yA8F?Pgg#y$1%P~A|5cqV*N%ME)?v?$QpZF9_w
zY$>L|{w$^qF%;Bkw7OQ28)+ueY#Ro^8wNDnU=5%fY(rvAjft?Pg$72OjkFD>&^8T;
zw1&jmxNz_gonYq02VztO$png#Q4;zYt(eH_n^DaX=QcC8Wx#w>(d)^eL?z#+Sj~99_e+
z1KaW)V{AB#2>}hjt6;;ZYGP^S-Pa_)xc^{K6sGgrI%_c_g>crl4T^du5%MSU?;50cPUYB|t-QFf#a`9uy(qIcH!5AEQ&MeWjrY5&?$$
zT0UoyLI=EejDB!>;U3x$-#|(b0^%_o-DbA+tEDkDqj&Vr6aqMlxV#(|o*lXb11DYz
zZwX|P8K1Y!26?4s=GP_PDBP2N3KQ;f!-4+c4f)uw)^xrGP{3+UxJHsYvU2AHf8Ap~Z`!A@rgJ&9rLn#}I6hJC9Krr^ycqg~<51H8uzy#DNY5l*H;QGCY)8)pY^PP%)JN**n`7O`Cm
zP!U*FLaDOUySW=I4-Jz&>X$A`-SorPGFzLIBKv`5X@d}iBIDX{&Q7uWCpf(nc2Ld<
zUis6?Az(qV(~$jwm_#09yg=ZGA|SNr$ZXh7ga|%yIyiIVS2vKFn8nF)?ehX
zJa6;f2Hkl7a%bs}ih%quAL5ApH
zD)Yb1*qY~Z2d2LH*~vj4AE{)`Fs(Am3hI*|ri+;C)kbkH+u@fr0n1^zguoah{9JuJ
z?;HMY-n@v?C)vrbc<0~|ji#Ci$iX%b(jbPMXjIf_8ygaBB-q*;Nj3o3lVIC>J3@^j
znHoZjAlnjc1ZgI~jR}lQXa>L=LI}Y|*xP7BX$^_70&Ej$A-2RKIMPF4k~42p(LL)4
z0_}9e6GBK`hgAfYiHbM_J+cz?4)Es6L6HBLEO5Hw)VBQ`@@XT4*Y><1hpj1(X~{Y8~9(wc3WfMqbWlXwA%^Sq#sb*{Q<@D#JgL%PZ}z>WxCk
zWgXAYt&Di5DsRlGDl!5oE}qO&tXst5}H^B!IfKacKQirFuXVlEm{L6E{8_!WC#N=hn-0!Vp+T@o}dBnsd2S5G)&A9y`fp@JbL|}DrHDyFD(L@MQ
z*vwnTi*`KZnc}$~-6l+PZx54H0Z+3Wb{6NQ+m87?!tg6|a6z2!Fk4}U161KPWBOZL
z;EkuzLFAobu`CDtgaCIIpkiJ4m7BQlT5qU#Mzt$>j5cd+j_J_%xm@QVAxoJ^VUs`a
zYEJm@JfKDu;&lQ94x-klW=gh=mMWIV%vzpMdB*MH9+ygn)n_`pm$U
zn;n>->ze?aBV)#h$^0l
zZk$~ETt;k4aq{I0L{^9kg)4Po5MHMHj$WSmEeO-JZJTv?bKJLU@rLx7uMJRditZaO
zSx3)*2t;?{*t`TdS;jCKC1
zoco?}k&{#-fr{oB;)LNkqNH-TSmV_w+IS;~$T0^OiEC+C25hKLa5RDOn9$wVV;Bh3Ol;wi4rdOyTL6lB=QonNr6k2HTr
zyPBi8@%}PuV03t5$g@KPDI?-0+^pw&4^h(fH8xeNwjuk`m4@!pZ5t~@O$2>w?smLA
zCcwPQ%#t9C5*r56Nw7fE1~l4iY#S3}X%VCfG?QQ?+5xc%um;m;O^LAyut*ycQGtm*
z4X>OZ~5POl<6iXhBi?N=L?{(4dB;_8r<-t5L~2|U}CgEZ|yk5t%tVU
zS)@F8dIA&k7J#n|LM0RkK>T2?<5N4!pOA46u}z!Tq7%FG*x}tpW7j+P@{Pvu%ax+;
z;~FB+&+tA)dYTs5LNr|1i9g}_2n2#?ToUev>u%E4vVbiZ93c=Es&+Nw-BfqY=(go<
z#{3}8jZk4?7F=L_FGdiA*K-Fsqx>tx1
zTj`!ocCjE4WYUWRo=`0Magflk?DI?j;3%SqB`Febl3?GhkSIY>bQ~^4K9M3IEwpA-
zz7;INo)hBWZ(r1BylxcAtRURa`)c=z=}P$6D{m&Dn{?`~Kzna~fA1P#>6&^-NA=T3
zlMoAvqkZw?y)17%57tVHa(7i}s
z6h8BOg-;H#2tqV`*Ky=RI0=s7As=JlD1ZjLW?3SO+x#_Nv#8_4kCS|Sec4i>Mxe==
z?xY0A93n9$4wqZeQ!Hr-Af?kKX
zAu@nKK@Z|KeU?Jw3T2o=^K67^v-g6wLT6`V8pds7WEe~kgrgZFRq8>#63oHlV;mt}
zr0WguYZwX;9K=iSOdG9={d@sUJsGbpvdCM#dIVS@06k&Va|&vY~UnJOSaQKKX1j4RubOlp193VDB{7H9L&;8%~PFr6-XBK~^og&d2BJ|qdW
zCkaAM?aa}bLW6r_T7W(|tb2yf_QS8UVV`K9>DcLiFLV^?qRH2J1qzU*C7i_9+VZbE
zTfN@dw-jq2v>^^%aP(Yy~dwzNo@YRM#a?~8~|7bGIN?Szy>NAUt(`ROz-hfv{`_(h!?KkOOJ5F$8E5
zNj8F-LP-&TCfLCPV1(F#v@G?27c`Jl5QK$@wkKLWzqIx~4109#j1Gv$W&8IzE=U|+
z9|!GR%k)qX0#yt`PL43N4fnS|)5+~jr7vxiq8Efo?qKPN0AUdL&>-NLtc3FGR)ik0
zvzPR1;$%;b%Xd0rbUy$$-~u2Yw~KX%0T4XG{R{97F@Z$iFG4~(m)*0O45qbWtLI<%
z%)%MlWe}`Bz}^)^s|hcgURkw5d4MWn(VNH`j^rv=;LfIREh-r!jIo
zq73yse9wMw-s&K7kQ8>d93kA;-ZRBq^p;_}MM2P-H8bu|-RFJ?PUTx?
z{ud9I#nRdt=H^%(C5{#N4-5Rop0FvSwaojff(@f1ou>tp9(++A!NgMpNza&Z0LD}=
z&Q7&1zzWrii;LeZq}z$a3zb>>c-2ago(p^~NAwx948p5klD2}-YlWXs9B+LiE7l_e
zE4EUF(3gkM{z6!f2K;nBIrRsJKrxK^6k+gBt~0$@b5GoC&P{nIf!{xLMlJ~hXvnU$
z3X$9X{&YV~98o9+ouH}d%v%%
zY+MrdMm9A3HzXkA6Ov438}WYL{bGxDjq#=&ja5sbBtFC3vzd&ks
zwGY(=hl)XUqZ;{^7o&NHs`mT&Jp(D`>KwUbGG}tmdc-n%pp~^!qErPGRUn~F!Hzb_
z;*nYu>A7cV3x5x52Snw}Nuha}`Lh)-bdF$S4qCpiGonOP6Ja_oq
zbyqUHv$8LbNQnt{y>y2@ZLcXh|Um(e`{?Xj|J~CRIlY7QtTh!kLWu
z(H?xL56kBxr&OC`iN-fp44G~|;?d{^S_@V(>r}il(@%{3^3Z;fEz%g
zNfD$d(chk7?Xl{%En~UY4*}-;6(7@YR$=gS!)Wnc#vj`+W{e{-
z!g1%!Rn;gR=-Y-(Bur*`vTUEISm8;G9v`v!^w=_X1_zW7n@Uk1&~3v8pdkfyCISm{
z=P6u@{$J59X$o|LTd{tnXt#CHfgh5RH|K?vYayW)Y#?j6l&LrCsCd&b(`OB99&Yj4
zx`9a;95QIh7@6j$P~8%yXmq6>`dF?B!GukEvidsF*6hKpYmW!Lih7bQfw1q)Tr@nd
z+Vf9)xTwGJ9wDdD)Ki)XSlt8t+fuREb=kGSEOfPVEo_FWp(Iamz~nchN#f7SaDM(_
zTd|G*a2MqEmcg9)DM`G9j2MRWC06jbr#6^C7+lMJuY;!Se41$6AX#$*~T)^4q
zFUo-Z`j{FVk;BArAZi|x>5BDwhiid&QL}QfWv&>kQm)!BHSaiI*1mZ~5CHJYg`8Q(
z`)TxYofGo%z4iRGSR|H{YbU}ZO9soqj)ECQZku=WjEqI)ljF9uW<2c+VN?(@$_b_9^{?_K
zZL$OpCQdgC*|ZwAq7uOjTfY*(pt#>(L`UJrM3`wix23|ns)<6
z#8+TO@XM4~J$OW|`Hhf29NPwuUiZ*zf!3OgX{k6JSCeX45QFK|->Hnr%DJGzhb~}1
zN)U&b>r>hCjXv1D#s!2ckn)6(`XuvH!3mGh6zpSQHoCMt(W@RHvNyv-Mgy^#G8e{}
zdLOpzaHhd4h@zFge}}|U*}QVhqqcO6>AiKR_|-Ew_N`}v<&5k~gT!Iq%75Q2=7)om
zvlEf;l7$7LieRYVbbrXwGJZpWaO~NPsUis32MA?%PXb{t>>LN!lIPPdH{7z5Dzd{m_N>Y%aB?wSU$z*6v!KD}wWN9YBup3Dx(AXidHX*dwHi0&hZ6?7X
zu_n?8z!L8M-bBOLLM+4^5sFC=Nwg->Y@yW(s<5IX5>IRew|~_H3xkXLkzHIWxCFe=
zlp$Jpa$!%{@5gG&OSqa4fMtr*Fk5b9;G5?UpN>GL3yLm^N@Jur&PhV9@c{WETwohD
zA!`27LYQFWA92N_KyeEJ=t^X>_WB{nRbu;oje!Jt=z=c5spdNQ?=Un7L{atL=W4sB
z1}U>X>BXk3QxJHi*r5o$C4IxEoP8J_Oo^vAw1Q0SoxrdmTt?fyW;KC)Sst1Q0|5~4
zCpq=-R#tHnU}ZFS$jQ)+?g3~^M6REaI7NoHs#dkSBT1G(7?hv@{2&CNY$I8iJUINV
z!S^_0?Sqzdd60l7B1Vo{AE!$!cOdu>gxZw&z&N|JUBJ0=n)ZE{F8VNCrMxUb4Q#xi
zVE2++M%x2*a&kV9V@wIN)tPW~%h5EWt=P6t+k>`a$MFiH7Ki7?gF@ILTSE{jU`gt4
z%Xof&uestp|7N+uKoO}(2u2k*8`!>5w_$u~z5Ht;@*Karp5b)nan79P#&!lOwIiZn
z+h~7PLZeNT?=?(zB9M~O(nn&3pVf6*Ro7&t{eu{O{r{yr{7L73Y#iZ<5`Og%fmG83
zsbx7ELAQ)y|2ff*MLt!xtg+dv#x{K
z+(Ar^ieoTpJ&{$v0Bk)zc|tks^;`EkP}hifOJQFbXg9K!O#ZNFZO2+~pHgoJ4uu4-
zbB6MaaME%zZXkw|lFDn8N1AS^RsS>6Gake3}39v(BV52}Hn39bKuXI3M(_kjS3A4`Hjr*~t`;T%L
z1w3d3n9i9vLRa|Aq$lx8n)C++BOP^265@37(umv7fp`wT
zDI=;e%ig#&2LCaT1eEDJQNB(VeX&v9>5BcSjCPGq=-I^r1bNDMr;izuu@E1lS+?X8
z%?$29@57h#XlRe;Fw^?*bx*YNzZ|_{%59R52NgsifoSu}=*Br&RJ;@uSl
z;)W)oq|(RjRL5qc1`*EjsK&(qIE>+jMbF|o{O6vpe-qyuJ(517XX7=K!{jwSQ3D3;
z#C_hNdW_vCMC2pdNGKs3@*mL?I6grxN^(dw@ON?mZu6x6=PnUjD1X4*_gd>U5H~fB
z_w_O-varblByn4kowdyrlJQ2dj3Z6aE$4cQH78qI=OvtZnqiP6j5cQUJh>{mk$W
zdsY}5ws+&exU?TrOa%zG&fb%gu^bHCSi0hV1NygANZ-~)3uL)vz
zyqi}r$H%E}28wp*?-vCWU+KY*e~r}?#VJP}IwOOUil|3oqiWLFj{j56GoIw|j%On}
zb1SBmV)}ZEo@LR?SWYNLHaI6N$++y!RI7pE43(Ns?XKY&%Jh@2J0YEGbqu9lpU*16E>BFmiN`fk|jkvI^jeHhr%!=U$tf*
zAZoH%M=f)3ArDcLvVz~69(y#jB!2B(O%q~%tD^lqKFK$WN-PnlpOKrT4Nw*S^1+_lz0n8e1eYBt8m$lFM4N(mGMi4ugQMXfTvtsxL<
z$^~4#j~wpt-W{kZ@L>q9`P!mln9sLkBqLif~3b>I9mWnH5B)op{J4zFT?~mvO3a=
zrmlvBxFDb|q-R)qeAz2>@^t&532MBtU?D-H2ijP;@W`RuJNtL$2SLYP{*M^D-^0T^
zK2We)kJv>T0~Be6p^7{D{SkkuI6H%U)S}<%hBc4(t%w>Lh@HXxo(lb4$uhKF+
z`^7lFka5k3ix9I#*$@*_R4|;nHqcQLmG=o^IDR+V*lB~_s}?ZiijI~1^4+DQ(l=;q
zqwAp(aM(9yQ*$-LF;Pryt1?HCGj%wRx$mr)*B3655L`r|1R(W{A`#{{dLrnT7@m5%
z2H|WL)Uy$NZczrYhMigb$IkC{4Jw|VHj+C9aJHPM!WJO(|LyamlYSkfE8N4^Y*jKVOd
z0FW`LxSc^7*@p=2o(QAGw4h)6g5Wi;wG#@ccs5+@^V$Vn6{8>}UJ%fR#m2NCj7SE=Ceji(b|ee0dy{@d$lGV6;Kg~7rrQwO8yMS3
zw&R0iY)z$rbzEoH%J;m+ILN==WF?-wm-zI4Z(wSaE4BmrtYLjMXaNIo^fu7c0SHF;
zs4z3pl4bgYyV8GxNBK)lugiWaXDWg-!F+y$qr>)~V+~7F_>>4fq`{CgTebBR;wW(Z
z(Mw%VNrT`3LK4s?s0^_HE6)*?3f}PQIeh9;5QL!U12+M#9s_K)Izce3C5Jc}fI&(Z
zXt2V--hETFx9-=X)K>&KzE!NI>u0=DjI$okTKSR(-dCp81~e27tp2~~8+W>S-uB%D
z@9Iv$UluN8?Y8|4|2%C3qPr!UFtXdD3@jQnRB=M=(YEJ-l3J{vacq!7Xg9;&4|E6w
zT2ANixdMq}WpR~j=hlS(KdRk1PO4%ZV-$M;sPl>OZeA}(gr92kT`r{CdqSE}4gL5N
zelGCs{UajnL%5NIgxGfGGDIWy!SP^@*y6L!3s<_mb01))kAD6L04IO(G!VPcBH%nR
z1uJ4Moo!jHvrMv;(4jQ4XJB3s1Qy0xbZpBXS*+*gQ0#!7n!sAay^)I#N3xWNGX#zx
zsC<6(-0$t|=R-Lg>W7?<$=wX&7F@{X9^oe_QX}`!v~TL&o4NARib^4}%nG2P`+X
z#Y{lV8>uJ(&M(R)(9}v&lIv7Y?K34%ufk<8ry8D!L^qokiTTb!X4-hV0S>x6La1fP
z2N|)Nqd%`mv5oI&DrKUqd`Opk*rv$blnLGSv(nLxa}`V-Kx$I_-w8Nj0mu^6P<;`F
z43WVU7~lv4jfg$l9pVR8E8?t(5yLQy*P=(N#Pke??;TG!eI@gmf3c6P7~m~*o4{=OCj>$BjI4Wk)`gJ2oKv!Z|-L~
zqQewoLmeR7V`!1Mh~sj8465XrYP6|q5>l3#UtxbB*A0W1%+B8SU!Z^XZKLx?`Z8@e
z;rD}}@J~w;y}`$Zun@i^Zjr_j^hN-XaXxm(r{cO#cp$?ru|qe;WW9=oG(gbgs_XxpR$2X4jNM`~if#%Moqxf#|`
z8*I1uDaT^iiQQ7cTYG&aFT|hZVQVJs^JCu4<_30cOOGoQ-6v&zKkQRi
zXSY>Q6(JO|!nh#<+)d1`&kL!h}JUZOJd5T==zcN&pwM`5U>LgGu}_OeQ<~s0QJK3R|{nceyOXeVYMko%Lvyx%UyDDOSq9Sn;+FvU8$yO-
z0!GwF7RfTiF#i4_zqUhqUuk1Ie9k~YXSW#(2Y3p2w|H*i7M#62t;qlZn=h^2%0BhC
zW?Qq_!Rs0^x^$E;RM_eg!YCQoMQ-DjEQVMRiPEjB#yPS^E_|+7XU{n1vJCNzSY}HK
zx1nSLCP)*Ys0f2RF%|ukrtORKi3(6n(O4r($;R1!i+0Gx8a^ogXRcy;Bjz1C
z-ZaP`wMw<5b0WL#Pb~BgE()%LgJs}>(eThDWBI?-UoB%_34dU~chdTw20PXsMBY+v
zQzHpqC{&`brY}K#T4+wzUav?Ny(U6Uwm%9lk-B1fwL6~^p?)ZJHlT_b1SoQg!EjZ`
zbT0?aV*W6j`YYcsko*Jln=^WBo9dh(%65OT2w)K&=k3r2FzSJ`5ZeY#Uue@=gN|Gw
z8f|Y%Ak0M$Qw{1i-l8Y{>?7%fNlekyTCnG8xmPS!B>vpDPymYJYhqdG)OUS&OGDw3
zHLOnDW>lbJD`Q;;7kJZ9n-OW2WG(OCYduT2M%}#Z$#VptVa89I0~3(8kPLt=vSbpv
zF9MD6PH#m}zXC*wNGHiqd5F*h(Awx&v=5lY`l4yK7X#pSZ>E;vQ{0UK4aasMuRv_8
z$Bt4q@c*bSz0APlUOuj5cMmd7Cxrg^DF^(#M*Uth!rr(|h-LE~vK>}khBag3X|~1A
zG87@_Pe6Fg$O)sXd^)c(8@k|B0uoe$_9c&*EpQ$J!kn!*bdeon+%l23`9cy(pZK)4
zud4Q%!Lc63^@=n-kKj10NsYrFIqe?ShGvG`4hYe?3k+hsIFaOi-hG3T-}f?LYPk^a
z-K0%D3h&Iu$_2iE*Npb%_oT*_i`L}4=ljyGGd8|f-OImk(wDwC+JB(PxP4KB$PdT)
z*GzxYjcEM#!nh#_ar!J!maP6P?7^H{|BU~th4c_Lw0!lq7^Q6{CevlVJdSTjLMsG0
zdR2P5p4D{S1vekFPpa>1nraO7yH~DXtlV_Rb-k2Lnyf0+dC4PlN@5@n#h7R7PaCZo
z@m^;!cWWu}T>-?*V=A79ziBMo?MYK~X2#?nr*G9%{wK4E{`Np~I3K_DvH)Rj+`As&inZ`D*_H_I^&`FSYO|do|8s0yGm_}Pp|ak
z$0BuAs;8E_`nU*E;=ciWDRD?8DLE=8jWE6nlyggXw~|B(L`nO*;#kl=W53|bI}6jFm&iGOSA51ruFX!)#Z!u(5=F)KjokG>TB2;f2k
zs&k4jAg!Q7hl5pqX2xr##|EArquomffq#D(_GBbEC)}4CDj=iwOS>3LWe-m!^j4`K
z%)GHlQ_M%x%r_*j6RNuLRXRmr{(P~FVW
zjU2()216WF1tH!KV+3bVZUT^K4Bmao-0DN$#|KMAhArSCe(rv;@0Zfo3NueV(c~2>
z+I0i3E&dX2y^`nXK&gp$)>4&2A%M$^r>4Z5Z@#u@+ltzY&zRl^(^3R|3tc3NEF$mC
z>cmf1W9-(3e+|vk@PI@@@U&@j1wh4KO?e13zLF|7AvfhE0|nyX3aCwznBg7h{9JVb
zDv><)(1i>w=Mmh>dC-O7f3x4XeMOns;i(|hnh~b7gp_Deu3Xp6)W>2NjmNey4%aU{
zbE`o*04hqDfAWPxzbf7x1hpM*t?{=Q(GVu4$)N8$`iGLLrIw~O02|!zf0DVhwT`ziXLq=Q#
z@%L0!X=kCn3~3U>T6J7SUXp}utA^U!w9hG|!PA-flhvL2Zl-V#d~PZ0x|A-5Cu2X)
z9agy43rV`!@ykh(IwRSsk1Isclq5;Ppbfl&I5dv+xf!8^(72`^cp6?Q#QGxa$crKR^WuC>NRd#K3>3ii;h
zHCRG1`SRgGCIsW`m;-Ofp`bMPW8YuO{c%^a9CdL@5RE!u;}i+5MR2LlYgp|hY=?=`5VzPk|+aC#r9l=(BCq=%aFp8^L0z)(PI3Vj)6u>)`Qo%m#HhkYL%0ufVdSZ?f{
zV_SgdV|$
zy5KU!ymn2=@FOU0Kdy$zf6uoFhR=%}uW`gUA^#S^Ln6R&U4)r)JCo$T<39DongJvLaE`+o8#A11X__GD>lPz
zGalWP5PVEp!}XAcFGAS*iQzmhM)))8)qn|tQb94(6L=*sxTDUsikEM4Y#y9LV_@!VmO?2hxqNz+8h#y4pOJs;309J)H
z_W~*DrZtM2p^pcVKhXF)EZiBD#=0iwD&2Fwga<&lwFDpmvYJ-1e6ajIVnaw(LXSUN
z5FFwSlFg{^PWKWEkWU6d&(W#8dor!`j2RTF0Uw>eyJ`M5ti{ctVX?I8MHinw^370gF(^sUBwz$Q
zARD5|J-ub*e}QEKmqmC*IdN^aO9ODQ5=WNoB&hU7%a&5}xiSQi7F`XP@L087c$#N(
zQ!Br_ePL;gTK)R0|FPftL=uN?Ckt&SX$-P5D(Cj{lRC`y!Z?j|
z`n)LHjGIGwPw-9{E={=0=Kpjx=NteJB(-hKOM{z#;Ck7uY@m&u24b4#hgk>_x8~8F
ze%h@jG=9?+3P!QV==AjS`nM+jB?(GWs~^p&LlUM!hb$H6AYa5R^VnHrB`83Ir2r);5CU@4a!*2Af8Vn9_F1le
zrE)G-RUIn47!nV2yUR#WDT8u4j@T(NsNp}K#8ODA*EaD688ODcY)}wW$<0!^gqyU{
zvE?Q-QS^SGbi>fcR|b3E=|Dmc2iWQh1V7Mq@k>Al?R#C_NK+lJU899z=uw}QTtCiS
z+K-HU@ON7Ji1aYAhMrEGYp{^=vET8P`q-f-da3CO!ALF!ng1}96_}vlVi!U%Ds#g(
zd~6+*rSv1xwaw6*ZO`E4`{p~Fw``A-V=TV;gM
zVQ;gxVzW0^k1eWcEj7o{?U1Gl%s5Qqt)Az%gNt`O2D(OaO^;&z+wPgAZTlw_cJ3Qb
zg+y9-{{v*CPWKOSkU^*w65mi~0f(s)69z{-;Q=T`Je}hg{1Z2m=;neBu=J(8_>`71
zaHn?Yc(MZvu+re#%{;es`{r!=UGrw0ON>XD>rOBOD9KD+bFOmAAD5Hm83AJ-^=j8v;C+!m72JuN_nZ@SUJI~oF1tuK_JFd=fSN?3&0m%YB_sr%c
zgVVghq!xwWZG4MT?81ma(LC*AkJ*4=gYrJ8fews#270h!rH?)-nSY|@o{GKXVqLM8
z%l}mGHz61^=2!@u4)b2qalCU&d)G7lyqd7IHa3M2qBOVoi;