Skip to content

GitLab

Commit afbe02ba authored by khalid's avatar khalid
Browse files

Rename params for RadSex tools

parent 22ae3904
Hide whitespace changes
Inline Side-by-side
tools/deseq2/deseq2.script.R
View file @ afbe02ba
......@@ -8,17 +8,20 @@ library(GenomicFeatures)
library(ggplot2)
samples = read.table(snakemake@config[["popmap_file"]],stringsAsFactors=F)$V1
conditions = read.table(snakemake@config[["popmap_file"]],stringsAsFactors=F)$V2
popmapData = read.table(snakemake@config[["popmap_file"]],stringsAsFactors=F)
colnames(popmapData)=c("samples","condition")
samples = popmapData$samples
conditions = popmapData$condition
#samples = c("DRR016125","DRR016126","DRR016127","DRR016128","DRR016129","DRR016130","DRR016131","DRR016132")
#files = file.path("/home/mbb/Documents/results2/kallisto/quant",samples,"abundance.h5")
files = sort(snakemake@input[["counts"]])
names(files) = samples
if (snakemake@config[["quantification"]] == "htseq_count")
if ("quantification" %in% names(snakemake@config) )
{
if(snakemake@config[["quantification"]] == "htseq_count")
{
print("Quantification was done with htseq-count !")
files = sort(snakemake@input[["counts"]])
names(files) = samples
sampleTable <- data.frame(
sampleName = samples,
fileName = files,
......@@ -30,28 +33,36 @@ if (snakemake@config[["quantification"]] == "htseq_count")
directory = "",
design = ~ condition
)
} else{ # counts where generated with one of these software "salmon", "sailfish", "alevin", "kallisto", "rsem", "stringtie"
if (snakemake@params[["tx2gene"]] == TRUE){
# conversion transcrits vers gènes
TxDb <- makeTxDbFromGFF(file = snakemake@params[["annotations"]])
k <- keys(TxDb, keytype = "TXNAME")
tx2gene <- select(TxDb, k, "GENEID", "TXNAME")
txi <- tximport(files, type = snakemake@config[["quantification"]], tx2gene = tx2gene)
} else {
txi <- tximport(files, type = snakemake@config[["quantification"]], txOut = TRUE)
}
groups = as.factor(conditions)
designMat <- model.matrix(~groups)
colnames(designMat)[2] = "condition"
sampleTable <- data.frame(condition = factor(conditions))
rownames(sampleTable) <- colnames(txi$counts)
dds <- DESeqDataSetFromTximport(txi, sampleTable, designMat)
} else{ # counts where generated with one of these software "salmon", "sailfish", "alevin", "kallisto", "rsem", "stringtie"
if (snakemake@params[["tx2gene"]] == TRUE){
# conversion transcrits vers gènes
TxDb <- makeTxDbFromGFF(file = snakemake@params[["annotations"]])
k <- keys(TxDb, keytype = "TXNAME")
tx2gene <- select(TxDb, k, "GENEID", "TXNAME")
txi <- tximport(files, type = snakemake@config[["quantification"]], tx2gene = tx2gene)
} else {
txi <- tximport(files, type = snakemake@config[["quantification"]], txOut = TRUE)
}
groups = as.factor(conditions)
designMat <- model.matrix(~groups)
colnames(designMat)[2] = "condition"
sampleTable <- data.frame(condition = factor(conditions))
rownames(sampleTable) <- colnames(txi$counts)
dds <- DESeqDataSetFromTximport(txi, sampleTable, designMat)
}
} else {
if("tsv_file" %in% names(snakemake@config)) #it's from a tsv file containing counts for samples
{
print("Counts are given in a file !")
countMatrix <- read.table (snakemake@config[["tsv_file"]], header=TRUE, sep="\t")
rownames(countMatrix) = countMatrix[,1]
countMatrix = countMatrix[,-1]
dds = DESeqDataSetFromMatrix(countData = countMatrix, colData = popmapData, design = ~ condition)
}
}
keep <- rowSums(counts(dds)) >= 10
dds <- dds[keep,]
......@@ -128,6 +139,8 @@ dev.off()
speciesData = snakemake@params[["biomaRtDataset"]]
#what if we have tarnscript ID or not ENSEMBL ID instead of gene ID ?????????
if (speciesData != "None") {
# Add annotation with biomaRt
library(biomaRt)
......@@ -137,7 +150,7 @@ if (speciesData != "None") {
species_mart <- useMart(biomart = "ENSEMBL_MART_ENSEMBL", host = "uswest.ensembl.org", dataset = speciesData)
#for the moment the top genes change to significant genes ?
#what if we have tarnscript ID instead of gene ID ?
ann <- getBM(filter="ensembl_gene_id",value=rownames(topGenes),attributes=c("ensembl_gene_id","description","transcript_length"),mart=species_mart)
go_ann <- getBM(filter="ensembl_gene_id",value=rownames(topGenes),attributes=c("ensembl_gene_id","description","go_id","name_1006","namespace_1003"),mart=species_mart)
......
tools/edger/edger.script.R
View file @ afbe02ba
......@@ -7,23 +7,44 @@ library(rhdf5)
library(GenomicFeatures)
library(ggplot2)
samples = read.table(snakemake@config[["group_file"]],stringsAsFactors=F)$V1
conditions = read.table(snakemake@config[["group_file"]],stringsAsFactors=F)$V2
files = snakemake@input[["counts"]]
names(files) = samples
if (snakemake@params[["tx2gene"]] == TRUE){
TxDb <- makeTxDbFromGFF(file = snakemake@params[["annotations"]])
k <- keys(TxDb, keytype = "TXNAME")
tx2gene <- select(TxDb, k, "GENEID", "TXNAME")
txi <- tximport(files, type = snakemake@config[["quantification"]], tx2gene = tx2gene)
} else {
txi <- tximport(files, type = snakemake@config[["quantification"]], txOut = TRUE)
popmapData = read.table(snakemake@config[["popmap_file"]],stringsAsFactors=F)
colnames(popmapData)=c("samples","condition")
samples = popmapData$samples
conditions = popmapData$condition
if ("quantification" %in% names(snakemake@config) )
{
if(snakemake@config[["quantification"]] == "htseq_count")
{
} else { #one of tools handled by tximport
files = snakemake@input[["counts"]]
names(files) = samples
if (snakemake@params[["tx2gene"]] == TRUE){
TxDb <- makeTxDbFromGFF(file = snakemake@params[["annotations"]])
k <- keys(TxDb, keytype = "TXNAME")
tx2gene <- select(TxDb, k, "GENEID", "TXNAME")
txi <- tximport(files, type = snakemake@config[["quantification"]], tx2gene = tx2gene)
} else {
txi <- tximport(files, type = snakemake@config[["quantification"]], txOut = TRUE)
}
cts <- txi$counts
dgelist <- DGEList(cts,group=as.factor(conditions))
}
} else {
if("tsv_file" %in% names(snakemake@config)) #it's from a tsv file containing counts for samples
{
print("Counts are given in a file !")
countMatrix <- read.table (snakemake@config[["tsv_file"]], header=TRUE, sep="\t")
rownames(countMatrix) = countMatrix[,1]
countMatrix = countMatrix[,-1]
dgelist = dgelist(counts = countMatrix, group = conditions)
}
}
cts <- txi$counts
dgelist <- DGEList(cts,group=as.factor(conditions))
#Filtering
countsPerMillion <- cpm(dgelist)
countCheck <- countsPerMillion > 1
......
tools/edger/edger.yaml
View file @ afbe02ba
......@@ -18,7 +18,7 @@
inputs:
[
{ name: counts, type: "tsv" , expand: true}, #will expand over counts for samples in SAMPLES
{ name: group_file, type: "popmap", file: "", description: "Path to tsv file with samples conditions"},
{ name: popmap_file, type: "popmap", file: "", description: "Path to tsv file with samples conditions"},
],
outputs:
[
......
tools/radsex_distrib/radsex_distrib.rule.snakefile
View file @ afbe02ba
......@@ -9,8 +9,8 @@ rule <step_name>__radsex_distrib:
# config["<step_name>__radsex_threads"]
params:
command = config["<step_name>__radsex_distrib_command"],
min_depth = config["<step_name>__radsex_min_depth"],
signif_threshold = config["<step_name>__radsex_signi_threshold"],
min_depth = config["<step_name>__radsex_distrib_min_depth"],
signif_threshold = config["<step_name>__radsex_distrib_signi_threshold"],
shell:
"{params.command} "+
......
tools/radsex_distrib/radsex_distrib.yaml
View file @ afbe02ba
......@@ -26,7 +26,7 @@
options:
[
{
name: "radsex_min_depth",
name: "radsex_distrib_min_depth",
prefix: "--min-depth",
value: 5,
min: 1,
......@@ -36,7 +36,7 @@
type: "numeric",
},
{
name: "radsex_signi_threshold",
name: "radsex_distrib_signi_threshold",
prefix: "--signif-threshold",
value: 0.05,
min: 0,
......
tools/radsex_map/radsex_map.rule.snakefile
View file @ afbe02ba
......@@ -10,10 +10,10 @@ rule <step_name>__radsex_map:
# config["<step_name>__radsex_threads"]
params:
command = config["<step_name>__radsex_map_command"],
min_quality = config["<step_name>__radsex_min_quality"],
min_frequency = config["<step_name>__radsex_min_frequency"],
min_depth = config["<step_name>__radsex_min_depth"],
signif_threshold = config["<step_name>__radsex_signif_threshold"],
min_quality = config["<step_name>__radsex_map_min_quality"],
min_frequency = config["<step_name>__radsex_map_min_frequency"],
min_depth = config["<step_name>__radsex_map_min_depth"],
signif_threshold = config["<step_name>__radsex_map_signif_threshold"],
shell:
"{params.command} "+
"--markers-file {input.markers_table} "+
......
tools/radsex_map/radsex_map.yaml
View file @ afbe02ba
......@@ -34,14 +34,14 @@
label: "Path to the genome file",
},
{
name: chromosome_file,
name: radsex_map_chromosome_file,
prefix: "",
type: input_file,
value: "",
label: "Path to the chromosome file",
},
{
name: "radsex_min_quality",
name: "radsex_map_min_quality",
prefix: "--min-quality",
value: 20,
min: 1,
......@@ -51,7 +51,7 @@
type: "numeric",
},
{
name: "radsex_min_frequency",
name: "radsex_map_min_frequency",
prefix: "--min-frequency",
value: 0.1,
min: 0,
......@@ -61,7 +61,7 @@
type: "numeric",
},
{
name: "radsex_min_depth",
name: "radsex_map_min_depth",
prefix: "--min-depth",
value: 5,
min: 1,
......@@ -71,7 +71,7 @@
type: "numeric",
},
{
name: "radsex_signif_threshold",
name: "radsex_map_signif_threshold",
prefix: "--signif-threshold",
value: 0.05,
min: 0,
......
tools/radsex_signif/radsex_signif.rule.snakefile
View file @ afbe02ba
......@@ -10,7 +10,7 @@ rule <step_name>__radsex_signif:
# config["<step_name>__radsex_threads"]
params:
command = config["<step_name>__radsex_signif_command"],
min_depth = config["<step_name>__radsex_min_depth"],
min_depth = config["<step_name>__radsex_signif_min_depth"],
signif_threshold = config["<step_name>__radsex_signif_threshold"],
shell:
"{params.command} "+
......
tools/radsex_signif/radsex_signif.yaml
View file @ afbe02ba
......@@ -28,7 +28,7 @@
options:
[
{
name: "radsex_min_depth",
name: "radsex_signif_min_depth",
prefix: "--min-depth",
value: 5,
min: 1,
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment