Title: | Automated functions for comparing various omic data from cbioportal.org |
---|---|
Description: | This package contains functions that allow analysing and comparing omic data across various cancers/cancer subgroups easily. So far, it is compatible with RNA-seq, microRNA-seq, microarray and methylation datasets that are stored on cbioportal.org. |
Authors: | Arman Shahrisa [aut, cre, cph], Maryam Tahmasebi Birgani [aut] |
Maintainer: | Arman Shahrisa <[email protected]> |
License: | Artistic-2.0 |
Version: | 1.29.0 |
Built: | 2024-11-14 05:52:52 UTC |
Source: | https://github.com/bioc/cbaf |
This function calculates frequency percentage, frequency ratio, mean value and median value of samples greater than specific cutoff in the selected study / subgroups of the study. Furthermore, it can looks for the five genes that contain the highest values in each study / study subgroup. It uses the data generated by obtainOneStudy()/obtainMultipleStudies() function.
automatedStatistics(submissionName, obtainedDataType = "multiple studies", calculate = c("frequencyPercentage", "frequencyRatio", "meanValue"), topGenes = TRUE, cutoff=NULL, round=TRUE)
automatedStatistics(submissionName, obtainedDataType = "multiple studies", calculate = c("frequencyPercentage", "frequencyRatio", "meanValue"), topGenes = TRUE, cutoff=NULL, round=TRUE)
submissionName |
a character string containing name of interest. It is used for naming the process. |
obtainedDataType |
a character string that specifies the type of input
data produced by the previous function. Two options are availabe:
|
calculate |
a character vector that containes the statistical procedures
users prefer the function to compute. The complete results can be obtained
by |
topGenes |
a logical value that, if set as TRUE, causes the function to
create three data.frame that contain the five top genes for each cancer. To
get all the three data.frames, "frequencyPercentage", "meanValue" and
"MedianValue" must have been included for |
cutoff |
a number used to limit samples to those that are greather than
this number (cutoff). The default value for methylation data is |
round |
a logical value that, if set to be |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
A new section in the BiocFileCache object that was created by one of
the obtainOneStudy() or obtainMultipleStudies() functions. It contains a list
that contains some or all of the following statistical measurements for
every gene group, based on what user has chosen: Frequency.Percentage
, Top.Genes.of.Frequency.Percentage
, Frequency.Ratio
,
Mean.Value
, Top.Genes.of.Mean.Value
, Median
,
Top.Genes.of.Median
.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(3,4)) automatedStatistics("test", obtainedDataType = "single study", calculate = c("frequencyPercentage", "frequencyRatio"))
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(3,4)) automatedStatistics("test", obtainedDataType = "single study", calculate = c("frequencyPercentage", "frequencyRatio"))
This function checks all the cancer studies that are registered in 'cbioportal.org' to examine whether or not they contain RNA-Seq, microRNA-Seq, microarray(mRNA), microarray(miRNA) and methylation data.
availableData(excelFileName, oneOfEach = FALSE)
availableData(excelFileName, oneOfEach = FALSE)
excelFileName |
a character string that is required to name the output and, if requested, excel file. |
oneOfEach |
a character string that is used to alter the function's
behavior to select unique one cancer of each type that contains data for the
requested technique. The default value is |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
An excel file that contains all the cancer studies versus available data types
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
This function removes the created databases in the cbaf package directory. This helps users to obtain the fresh data from cbioportal.org.
cleanDatabase(databaseNames = NULL)
cleanDatabase(databaseNames = NULL)
databaseNames |
a character vector that contains name of databases that
will be removed. The default value in |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
prints the number of removed databases.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
cleanDatabase(databaseNames = "Whole")
cleanDatabase(databaseNames = "Whole")
This function can prepare heatmap for 'frequency percentage', 'mean value' and 'median value' data provided by automatedStatistics() function.
heatmapOutput(submissionName, shortenStudyNames = TRUE, geneLimit = 50, rankingMethod = "variation", heatmapFileFormat = "TIFF", resolution = 600, RowCex = "auto", ColCex = "auto", heatmapMargines = "auto", rowLabelsAngle = 0, columnLabelsAngle = 45, heatmapColor = "RdBu", reverseColor = TRUE, transposedHeatmap = FALSE, simplifyBy = FALSE, genesToDrop = FALSE)
heatmapOutput(submissionName, shortenStudyNames = TRUE, geneLimit = 50, rankingMethod = "variation", heatmapFileFormat = "TIFF", resolution = 600, RowCex = "auto", ColCex = "auto", heatmapMargines = "auto", rowLabelsAngle = 0, columnLabelsAngle = 45, heatmapColor = "RdBu", reverseColor = TRUE, transposedHeatmap = FALSE, simplifyBy = FALSE, genesToDrop = FALSE)
submissionName |
a character string containing name of interest. It is used for naming the process. |
shortenStudyNames |
a logical vector. If the value is set as TRUE, function will try to remove the last part of the cancer names aiming to shorten them. The removed segment usually contains the name of scientific group that has conducted the experiment. |
geneLimit |
if large number of genes exist in at least one gene group,
this option can be used to limit the number of genes that are shown on
heatmap. For instance, |
rankingMethod |
a character value that determines how genes will be
ranked prior to drawing heatmap. |
heatmapFileFormat |
This option enables the user to select the desired
image file format of the heatmaps. The default value is |
resolution |
a number. This option can be used to adjust the resolution of the output heatmaps as 'dot per inch'. The defalut value is 600. |
RowCex |
a number that specifies letter size in heatmap row names,
which ranges from 0 to 2. If |
ColCex |
a number that specifies letter size in heatmap column names,
which ranges from 0 to 2. If |
heatmapMargines |
a numeric vector that is used to set heatmap margins.
If |
rowLabelsAngle |
a number that determines the angle with which the gene names are shown in heatmaps. The default value is 0 degree. |
columnLabelsAngle |
a number that determines the angle with which the studies/study subgroups names are shown in heatmaps. The default value is 45 degree. |
heatmapColor |
a character string that defines heatmap color. The
default value is |
reverseColor |
a logical value that reverses the color gradiant for heatmap(s). |
transposedHeatmap |
a logical value that transposes heatmap rows to columns and vice versa. |
simplifyBy |
a number that tells the function to change the values
smaller than that to zero. The purpose behind this option is to facilitate
recognizing candidate genes. Therefore, it is not suited for publications. It
has the same unit as |
genesToDrop |
a character vector. Gene names within this vector will be
omitted from heatmap.The default value is |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
Based on preference, three heatmaps for "Frequency.Percentage"
, "Mean.Value"
and "Median.value"
can be generated. If more
than one group of genes are entered, output for each group will be strored in
a separate sub-directory.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(3,4)) automatedStatistics("test", obtainedDataType = "single study", calculate = c("frequencyPercentage", "frequencyRatio")) heatmapOutput(submissionName = "test")
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(3,4)) automatedStatistics("test", obtainedDataType = "single study", calculate = c("frequencyPercentage", "frequencyRatio")) heatmapOutput(submissionName = "test")
This function Obtains the requested data for the given genes across multiple cancer studies. It can check whether or not all genes are included in cancer studies and, if not, looks for the alternative gene names.
obtainMultipleStudies(genesList, submissionName, studiesNames, desiredTechnique, cancerCode = FALSE, validateGenes = TRUE)
obtainMultipleStudies(genesList, submissionName, studiesNames, desiredTechnique, cancerCode = FALSE, validateGenes = TRUE)
genesList |
a list that contains at least one gene group |
submissionName |
a character string containing name of interest. It is used for naming the process. |
studiesNames |
a character vector or a matrix that containes desired
cancer names. The character vector containes standard
names of cancer studies that can be found on cbioportal.org, such as
|
desiredTechnique |
a character string that is one of the following
techniques: |
cancerCode |
a logical value that tells the function to use cbioportal
abbreviated cancer names instead of complete cancer names, if set to be
|
validateGenes |
a logical value that, if set to be |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
a BiocFileCach object that contains the obtained data without further
processing. Name of the object is combination of bfc_
and
submissionName
. Inside it, there is a section for the obtained data,
which is stored as a list. At first level, this list is subdivided into
diferent groups based on the list of genes that user has given the function,
then each gene group itself contains one matrix for every cancer study.
Additonally, if validateGenes = TRUE
, another section that contains
gene validation results will be created in the BiocFileCach object.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) studies <- c("Acute Myeloid Leukemia (TCGA, Provisional)", "Adrenocortical Carcinoma (TCGA, Provisional)", "Bladder Urothelial Carcinoma (TCGA, Provisional)", "Brain Lower Grade Glioma (TCGA, Provisional)", "Breast Invasive Carcinoma (TCGA, Provisional)") obtainMultipleStudies(genes, "test2", studies, "RNA-Seq")
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) studies <- c("Acute Myeloid Leukemia (TCGA, Provisional)", "Adrenocortical Carcinoma (TCGA, Provisional)", "Bladder Urothelial Carcinoma (TCGA, Provisional)", "Brain Lower Grade Glioma (TCGA, Provisional)", "Breast Invasive Carcinoma (TCGA, Provisional)") obtainMultipleStudies(genes, "test2", studies, "RNA-Seq")
This function Obtains the requested data for the given genes across multiple subgroups of a cancer. It can check whether or not all genes are included in subgroups of a cancer study and, if not, looks for the alternative gene names.
obtainOneStudy(genesList, submissionName, studyName, desiredTechnique, desiredCaseList = FALSE, validateGenes = TRUE)
obtainOneStudy(genesList, submissionName, studyName, desiredTechnique, desiredCaseList = FALSE, validateGenes = TRUE)
genesList |
a list that contains at least one gene group |
submissionName |
a character string containing name of interest. It is used for naming the process. |
studyName |
a character string showing the desired cancer name. It is an
standard cancer study name that can be found on cbioportal.org, such as
|
desiredTechnique |
a character string that is one of the following
techniques: |
desiredCaseList |
a numeric vector that contains the index of desired
cancer subgroups, assuming the user knows index of desired subgroups. If not,
desiredCaseList is set to |
validateGenes |
a logical value that, if set to be 'TRUE', causes the function to check each cancer study to find whether or not each gene has a record. If a cancer doesn't have a record for specific gene, function looks for alternative gene names that cbioportal might use instead of the given gene name. |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
a BiocFileCach object that contains the obtained data without further processing. Name of the object is combination of 'bfc_' and submissionName. Inside it, there is a section for the obtained data, which is stored as a list. At first level, this list is subdivided into diferent groups based on the list of genes that user has given the function, then each gene group itself contains one matrix for every study subgroup. Additonally, if validateGenes = TRUE, another section that contains gene validation results will be created in the BiocFileCach object.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(2,3,4,5))
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(2,3,4,5))
This function Obtains the requested data for the given genes across multiple cancer studie. It can check whether or not all genes are included in cancer studies and and, if not, looks for the alternative gene names. Then it calculates frequency percentage, frequency ratio, mean value and median value of samples greather than specific value in the selected cancer studies. Furthermore, it looks for the five genes that comprise the highest values in each cancer study.
processMultipleStudies(genesList, submissionName, studiesNames, desiredTechnique, cancerCode = FALSE, validateGenes = TRUE, calculate = c("frequencyPercentage", "frequencyRatio", "meanValue"), cutoff=NULL, round=TRUE, topGenes = TRUE, shortenStudyNames = TRUE, geneLimit = 50, rankingMethod = "variation", heatmapFileFormat = "TIFF", resolution = 600, RowCex = "auto", ColCex = "auto", heatmapMargines = "auto", rowLabelsAngle = 0, columnLabelsAngle = 45, heatmapColor = "RdBu", reverseColor = TRUE, transposedHeatmap = FALSE, simplifyBy = FALSE, genesToDrop = FALSE, transposeResults = FALSE, downloadOnServer = FALSE)
processMultipleStudies(genesList, submissionName, studiesNames, desiredTechnique, cancerCode = FALSE, validateGenes = TRUE, calculate = c("frequencyPercentage", "frequencyRatio", "meanValue"), cutoff=NULL, round=TRUE, topGenes = TRUE, shortenStudyNames = TRUE, geneLimit = 50, rankingMethod = "variation", heatmapFileFormat = "TIFF", resolution = 600, RowCex = "auto", ColCex = "auto", heatmapMargines = "auto", rowLabelsAngle = 0, columnLabelsAngle = 45, heatmapColor = "RdBu", reverseColor = TRUE, transposedHeatmap = FALSE, simplifyBy = FALSE, genesToDrop = FALSE, transposeResults = FALSE, downloadOnServer = FALSE)
genesList |
a list that contains at least one gene group |
submissionName |
a character string containing name of interest. It is used for naming the process. |
studiesNames |
a character vector or a matrix that containes desired
cancer names. The character vector containes standard names of cancer studies
that can be found on cbioportal.org, such as
|
desiredTechnique |
a character string that is one of the following
techniques: |
cancerCode |
a logical value that tells the function to use cbioportal
abbreviated cancer names instead of complete cancer names, if set to be
|
validateGenes |
a logical value that, if set to be |
calculate |
a character vector that containes the statistical procedures
users prefer the function to compute. The complete results can be obtained
by |
cutoff |
a number used to limit samples to those that are greather than
this number (cutoff). The default value for methylation data is |
round |
a logical value that, if set to be |
topGenes |
a logical value that, if set as |
shortenStudyNames |
a logical vector. If the value is set as TRUE, function will try to remove the last part of the cancer names aiming to shorten them. The removed segment usually contains the name of scientific group that has conducted the experiment. |
geneLimit |
if large number of genes exist in at least one gene group,
this option can be used to limit the number of genes that are shown on
heatmap. For instance, |
rankingMethod |
a character value that determines how genes will be
ranked prior to drawing heatmap. |
heatmapFileFormat |
This option enables the user to select the desired
image file format of the heatmaps. The default value is |
resolution |
a number. This option can be used to adjust the resolution of the output heatmaps as 'dot per inch'. The defalut value is 600. |
RowCex |
a number that specifies letter size in heatmap row names,
which ranges from 0 to 2. If |
ColCex |
a number that specifies letter size in heatmap column names,
which ranges from 0 to 2. If |
heatmapMargines |
a numeric vector that is used to set heatmap margins.
If |
rowLabelsAngle |
a number that determines the angle with which the gene names are shown in heatmaps. The default value is 0 degree. |
columnLabelsAngle |
a number that determines the angle with which the studies/study subgroups names are shown in heatmaps. The default value is 45 degree. |
heatmapColor |
a character string that defines heatmap color. The
default value is |
reverseColor |
a logical value that reverses the color gradiant for heatmap(s). |
transposedHeatmap |
a logical value that transposes heatmap rows to columns and vice versa. |
simplifyBy |
a number that tells the function to change the values
smaller than that to zero. The purpose behind this option is to facilitate
recognizing candidate genes. Therefore, it is not suited for publications. It
has the same unit as |
genesToDrop |
a character vector. Gene names within this vector will be
omitted from heatmap.The default value is |
transposeResults |
a logical value that enables the function to replace the columns and rows of data. |
downloadOnServer |
a logical value that enables the function to download data on server and store the database as a zip file. This zip file can be imported later into R locally for further steps that can't run on server. |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
a BiocFileCache object that containes some or all of the following
groups, based on what user has chosen: obtainedData
,
validationResults
, frequencyPercentage
,
Top.Genes.of.Frequency.Percentage
, frequencyRatio
,
meanValue
, Top.Genes.of.Mean.Value
, medianValue
,
Top.Genes.of.Median.Value
. It also saves these results in one excel
files for convenience. Based on preference, three heatmaps for frequency
percentage, mean value and median can be generated. If more than one group of
genes is entered, output for each group will be strored in a separate
sub-directory.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) studies <- c("Acute Myeloid Leukemia (TCGA, Provisional)", "Adrenocortical Carcinoma (TCGA, Provisional)", "Bladder Urothelial Carcinoma (TCGA, Provisional)", "Brain Lower Grade Glioma (TCGA, Provisional)", "Breast Invasive Carcinoma (TCGA, Provisional)") processMultipleStudies(genes, "test2", studies, "RNA-Seq", calculate = c("frequencyPercentage", "frequencyRatio"), heatmapMargines = c(16,10), RowCex = 1, ColCex = 1)
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) studies <- c("Acute Myeloid Leukemia (TCGA, Provisional)", "Adrenocortical Carcinoma (TCGA, Provisional)", "Bladder Urothelial Carcinoma (TCGA, Provisional)", "Brain Lower Grade Glioma (TCGA, Provisional)", "Breast Invasive Carcinoma (TCGA, Provisional)") processMultipleStudies(genes, "test2", studies, "RNA-Seq", calculate = c("frequencyPercentage", "frequencyRatio"), heatmapMargines = c(16,10), RowCex = 1, ColCex = 1)
This function Obtains the requested data for the given genes across multiple subgroups of a cancer. It can check whether or not all genes are included in subgroups of a cancer study and, if not, looks for the alternative gene names. Then it calculates frequency percentage, frequency ratio, mean value and median value of samples greather than specific value in the selected subgroups of the cancer. Furthermore, it looks for the five genes that comprise the highest values in each cancer study subgroup.
processOneStudy(genesList, submissionName, studyName, desiredTechnique , desiredCaseList = FALSE, validateGenes = TRUE, calculate = c("frequencyPercentage", "frequencyRatio", "meanValue"), cutoff=NULL, round=TRUE, topGenes = TRUE, shortenStudyNames = TRUE, geneLimit = 50, rankingMethod = "variation", heatmapFileFormat = "TIFF", resolution = 600, RowCex = "auto", ColCex = "auto", heatmapMargines = "auto", rowLabelsAngle = 0, columnLabelsAngle = 45, heatmapColor = "RdBu", reverseColor = TRUE, transposedHeatmap = FALSE, simplifyBy = FALSE, genesToDrop = FALSE, transposeResults = FALSE)
processOneStudy(genesList, submissionName, studyName, desiredTechnique , desiredCaseList = FALSE, validateGenes = TRUE, calculate = c("frequencyPercentage", "frequencyRatio", "meanValue"), cutoff=NULL, round=TRUE, topGenes = TRUE, shortenStudyNames = TRUE, geneLimit = 50, rankingMethod = "variation", heatmapFileFormat = "TIFF", resolution = 600, RowCex = "auto", ColCex = "auto", heatmapMargines = "auto", rowLabelsAngle = 0, columnLabelsAngle = 45, heatmapColor = "RdBu", reverseColor = TRUE, transposedHeatmap = FALSE, simplifyBy = FALSE, genesToDrop = FALSE, transposeResults = FALSE)
genesList |
a list that contains at least one gene group |
submissionName |
a character string containing name of interest. It is used for naming the process. |
studyName |
a character string showing the desired cancer name. It is an
standard cancer study name that can be found on cbioportal.org, such as
|
desiredTechnique |
a character string that is one of the following
techniques: |
desiredCaseList |
a numeric vector that contains the index of desired
cancer subgroups, assuming the user knows index of desired subgroups. If not,
desiredCaseList is set to |
validateGenes |
a logical value that, if set to be |
calculate |
a character vector that containes the statistical procedures
users prefer the function to compute. The complete results can be obtained
by |
cutoff |
a number used to limit samples to those that are greather than
specific number (cutoff). The default value for methylation data is 0.8 while
gene expression studies use default value of 2. For methylation studies, it
is |
round |
a logical value that, if set to be |
topGenes |
a logical value that, if set as |
shortenStudyNames |
a logical vector. If the value is set as
|
geneLimit |
if large number of genes exist in at least one gene group,
this option can be used to limit the number of genes that are shown on
heatmap. For instance, |
rankingMethod |
a character value that determines how genes will be
ranked prior to drawing heatmap. |
heatmapFileFormat |
This option enables the user to select the desired
image file format of the heatmaps. The default value is |
resolution |
a number. This option can be used to adjust the resolution of the output heatmaps as 'dot per inch'. The defalut value is 600. |
RowCex |
a number that specifies letter size in heatmap row names,
which ranges from 0 to 2. If |
ColCex |
a number that specifies letter size in heatmap column names,
which ranges from 0 to 2. If |
heatmapMargines |
a numeric vector that is used to set heatmap margins.
If |
rowLabelsAngle |
a number that determines the angle with which the gene names are shown in heatmaps. The default value is 0 degree. |
columnLabelsAngle |
a number that determines the angle with which the studies/study subgroups names are shown in heatmaps. The default value is 45 degree. |
heatmapColor |
a character string that defines heatmap color. The
default value is |
reverseColor |
a logical value that reverses the color gradiant for heatmap(s). |
transposedHeatmap |
a logical value that transposes heatmap rows to columns and vice versa. |
simplifyBy |
a number that tells the function to change the values
smaller than that to zero. The purpose behind this option is to facilitate
recognizing candidate genes. Therefore, it is not suited for publications. It
has the same unit as |
genesToDrop |
a character vector. Gene names within this vector will be
omitted from heatmap.The default value is |
transposeResults |
a logical value that enables the function to replace the columns and rows of data. |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
a BiocFileCache object that containes some or all of the following
groups, based on what user has chosen: ObtainedData
,
validationResults
, frequencyPercentage
,
Top.Genes.of.Frequency.Percentage
, frequencyRatio
,
meanValue
, Top.Genes.of.Mean.Value
, medianValue
,
Top.Genes.of.Median.Value
. It also saves these results in one excel
files for convenience. Based on preference, three heatmaps for frequency
percentage, mean value and median can be generated. If more than one group of
genes is entered, output for each group will be strored in a separate
sub-directory.
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) processOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(2,3,4,5), calculate = c("frequencyPercentage", "frequencyRatio"), heatmapMargines = c(16, 10), RowCex = 1, ColCex = 1)
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) processOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(2,3,4,5), calculate = c("frequencyPercentage", "frequencyRatio"), heatmapMargines = c(16, 10), RowCex = 1, ColCex = 1)
This function generates excel files containing gene validation and all selected statistical methods. It uses outputs of obtainOneStudy()/obtainMultipleStudies() and automatedStatistics() functions.
xlsxOutput(submissionName, transposeResults = FALSE)
xlsxOutput(submissionName, transposeResults = FALSE)
submissionName |
a character string containing name of interest. It is used for naming the process. |
transposeResults |
a logical value that enables the function to replace the columns and rows of data. |
Package: | cbaf |
Type: | Package |
Version: | 1.27.4 |
Date: | 2024-06-30 |
License: | Artistic-2.0 |
It generates one excel file for each gene group. This excel file contains output of automatedStatistics() and validation result from output of either obtainOneStudy() or obtainMultipleStudies().
Arman Shahrisa, [email protected] [maintainer, copyright holder]
Maryam Tahmasebi Birgani, [email protected]
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(3,4)) automatedStatistics("test", obtainedDataType = "single study", calculate = c("frequencyPercentage", "frequencyRatio")) xlsxOutput("test")
genes <- list(K.demethylases = c("KDM1A", "KDM1B", "KDM2A", "KDM2B", "KDM3A", "KDM3B", "JMJD1C", "KDM4A"), K.methyltransferases = c("SUV39H1", "SUV39H2", "EHMT1", "EHMT2", "SETDB1", "SETDB2", "KMT2A", "KMT2A")) obtainOneStudy(genes, "test", "Breast Invasive Carcinoma (TCGA, Cell 2015)", "RNA-Seq", desiredCaseList = c(3,4)) automatedStatistics("test", obtainedDataType = "single study", calculate = c("frequencyPercentage", "frequencyRatio")) xlsxOutput("test")