Having made that file we can run DS (this assumes you have successfully installed the leafcutter R package as described under Installation above)
Having made that file we can run DS (this assumes you have successfully installed the leafcutter R package as described under Installation above)
../scripts/leafcutter_ds.R --num_threads 4 --exons_file=../leafcutter/data/gencode19_exons.txt.gz ../example_data/testYRIvsEU_perind_numers.counts.gz ../example_data/test_diff_intron.txtRunning ../scripts/leafcutter_ds.R -h will give usage info for this script. Here have included an exons_file: this is optional but allows LeafCutter to annotate clusters according to which gene they correspond to. An example exon file for hg19 derived from GENCODE v19 is included with LeafCutter. The exons file should have the following columns: chr, start, end, strand, gene_name. We provide a helper script gtf_to_exons.R which will convert a .gtf to the required format. I’m also hosting the exon file for GRCh38 GENCODE v.26 here.
Two tab-separated text files are output:
diff --git a/docs/articles/Visualization.html b/docs/articles/Visualization.html index 120b2c6..0453f50 100644 --- a/docs/articles/Visualization.html +++ b/docs/articles/Visualization.html @@ -112,6 +112,7 @@This step annotates each intron in each cluster at a given false discovery rate.
-Rscript prepare_results.R --iFolder <iFolder> \
- --oFolder <oFolder> \
- --support <support> \
- --annotation_code <annotation_code> \
- --code <code> \
- --FDR <FDR> \ iFolder - the folder that contains the results of the differential intron excision analysis. The script assumes all results share the same path and code.
-oFolder - the folder where the results will be written.
-support - the same support file as used for the previous analysis.
-annotation_code - as before.
-code - the same code used for the rest of the analysis, eg testYRIvsEU
FDR - the benjamini-hochberg false discovery rate with which to filter the results.
-This will create the Rdata object <oFolder>/results/<code>_results.Rdata. The file ‘prepare_example.sh’ shows how this would be done for the example dataset if you wanted to rebuild ‘Brain_vs_Heart_results.Rdata’.
This step annotates each intron in each cluster at a given false discovery rate and generates a single .RData file with everything the shiny app needs to run.
./prepare_results.R [options] <name>_perind_numers.counts.gz <name>_cluster_significance.txt <name>_effect_sizes.txt annotation_code ** leafcutter_ds.R. ** leafcutter_ds.R. ** leafcutter_ds.R. ** annotation_code ** will be something like annotation_codes/gencode_hg19/gencode_hg19 (see above)
Options:
+ -o OUTPUT, --output=OUTPUT
+ The output file that will be created ready for loading by run_leafviz.R [leafviz.RData]
+
+ -m META_DATA_FILE, --meta_data_file=META_DATA_FILE
+ The support file used in the differential splicing analysis. Columns should be file name and condition
+
+ -f FDR, --FDR=FDR
+ the adjusted p value threshold to use [0.05]
+
+ -c CODE, --code=CODE
+ A name for this analysis (will be available in leafviz through the Summary tab). [leafcutter_ds]
+
+ -h, --help
+ Show this help message and exitThis will create the Rdata file wherever --output is pointed. The file ‘prepare_example.sh’ shows how this would be done for the example dataset if you wanted to rebuild ‘Brain_vs_Heart_results.Rdata’.
As a concrete example, let’s assume you just ran the example at Usage and you’re currently in the /example_data subdirectory, then running
../leafviz/prepare_results.R --meta_data_file test_diff_intron.txt --code leafcutter testYRIvsEU_perind_numers.counts.gz leafcutter_ds_cluster_significance.txt leafcutter_ds_effect_sizes.txt ../leafviz/annotation_codes/gencode_hg19/gencode_hg19should create an leafviz.RData file.
cd leafvis/
-Rscript run_leafvis.R <oFolder>/results/<code>_results.RdataFrom the example_data directory.
cd ../leafviz/
+./run_leafviz.R ../example_data/leafviz.RdataThis will load in the Rdata object and run the LeafCutter Visualisation App in your browser.
Parse output of differential_splicing and make a per cluster results table
Parse output of differential_splicing and make a per cluster results table
cluster_results_table(results)@@ -100,7 +100,7 @@
From differential_splicing
Can be one of "smart", "random". smart uses an method of moments estimator to get a reasonable initialization. The seed for "random" can be set through the ... arguments passed to rstan::optimizing.
@param ... will be passed on the rstan::optimizing, so can be used for example to set the algorithm used (default is LBFGS) or the random seed if random initialization is requested.
Used to protect against colinear covariates.
will be passed on the rstan::optimizing, so can be used for example to set the algorithm used (default is LBFGS) or the random seed if random initialization is requested.
We recommend using dirichlet_multinomial_anova_mc instead.
We recommend using dirichlet_multinomial_anova_mc instead.
dirichlet_multinomial_glm(x, y, concShape = 1.0001, concRate = 1e-04)diff --git a/docs/reference/leaf_cutter_effect_sizes.html b/docs/reference/leaf_cutter_effect_sizes.html index b8174c7..0ba6a1c 100644 --- a/docs/reference/leaf_cutter_effect_sizes.html +++ b/docs/reference/leaf_cutter_effect_sizes.html @@ -100,7 +100,7 @@
From differential_splicing
\
- --FDR \
+./prepare_results.R [options] _perind_numers.counts.gz _cluster_significance.txt _effect_sizes.txt annotation_code
```
+** _perind_numers.counts.gz ** The same counts matrix you gave to `leafcutter_ds.R`.
+** _cluster_significance.txt ** The cluster significant table output by `leafcutter_ds.R`.
+** _effect_sizes.txt ** The per junction effect size table output by `leafcutter_ds.R`.
+** annotation_code ** will be something like `annotation_codes/gencode_hg19/gencode_hg19` (see above)
-**iFolder** - the folder that contains the results of the differential intron excision analysis. The script assumes all results share the same path and **code**.
-
-**oFolder** - the folder where the results will be written.
+Options:
+```
+ -o OUTPUT, --output=OUTPUT
+ The output file that will be created ready for loading by run_leafviz.R [leafviz.RData]
-**support** - the same support file as used for the previous analysis.
+ -m META_DATA_FILE, --meta_data_file=META_DATA_FILE
+ The support file used in the differential splicing analysis. Columns should be file name and condition
-**annotation_code** - as before.
+ -f FDR, --FDR=FDR
+ the adjusted p value threshold to use [0.05]
-**code** - the same code used for the rest of the analysis, eg `testYRIvsEU`
+ -c CODE, --code=CODE
+ A name for this analysis (will be available in leafviz through the Summary tab). [leafcutter_ds]
-**FDR** - the benjamini-hochberg false discovery rate with which to filter the results.
+ -h, --help
+ Show this help message and exit
+```
+This will create the Rdata file wherever `--output` is pointed. The file 'prepare_example.sh' shows how this would be done for the example dataset if you wanted to rebuild 'Brain_vs_Heart_results.Rdata'.
-This will create the Rdata object `/results/_results.Rdata`. The file 'prepare_example.sh' shows how this would be done for the example dataset if you wanted to rebuild 'Brain_vs_Heart_results.Rdata'.
+As a concrete example, let's assume you just ran the example at [Usage](./Usage.html) and you're currently in the `/example_data` subdirectory, then running
+```
+../leafviz/prepare_results.R --meta_data_file test_diff_intron.txt --code leafcutter testYRIvsEU_perind_numers.counts.gz leafcutter_ds_cluster_significance.txt leafcutter_ds_effect_sizes.txt ../leafviz/annotation_codes/gencode_hg19/gencode_hg19
+```
+should create an `leafviz.RData` file.
### Step 2. Visualise the results
+From the `example_data` directory.
```
-cd leafvis/
-Rscript run_leafvis.R /results/_results.Rdata
+cd ../leafviz/
+./run_leafviz.R ../example_data/leafviz.Rdata
```
-
-This will load in the Rdata object and run the LeafCutter Visualisation App in your browser.
+This will load in the Rdata object and run the LeafCutter Visualisation App in your browser.
### Features
diff --git a/leafviz/prepare_results.R b/leafviz/prepare_results.R
index 6fae0a9..230e9b5 100755
--- a/leafviz/prepare_results.R
+++ b/leafviz/prepare_results.R
@@ -5,21 +5,14 @@
## and prepare for visualisation
library(optparse)
-require(leafcutter)
-library(data.table)
-library(stringr)
-library(dplyr)
-library(magrittr)
-
-options(echo=TRUE)
option_parser=OptionParser(
- usage="%prog [options] _perind_numers.counts.gz _cluster_significance.txt _effect_sizes.txt annotation_code \nThe annotation_code should be something like annotation_codes/annotation_codes/gencode_hg19/gencode_hg19.\nresults.RData",
+ usage="%prog [options] _perind_numers.counts.gz _cluster_significance.txt _effect_sizes.txt annotation_code \nThe annotation_code should be something like annotation_codes/gencode_hg19/gencode_hg19.",
option_list=list(
- make_option( c("-o","--output"), default="leafviz.RData", help="The output file that will be created ready for loading by run_leafviz.R"),
- make_option( "--meta_data_file", default=NULL, help="The support file used in the differential splicing analysis. Columns should be file name and condition"),
- make_option( c("-f","--FDR"), default=0.05, help = "the adjusted p value threshold to use"),
- make_option( "--code", default="leafcutter_ds", help = "A name for this analysis (will be available in leafviz through the Summary tab."))
+ make_option( c("-o","--output"), default="leafviz.RData", help="The output file that will be created ready for loading by run_leafviz.R [%default]"),
+ make_option( c("-m","--meta_data_file"), default=NULL, help="The support file used in the differential splicing analysis. Columns should be file name and condition"),
+ make_option( c("-f","--FDR"), default=0.05, help = "the adjusted p value threshold to use [%default]"),
+ make_option( c("-c","--code"), default="leafcutter_ds", help = "A name for this analysis (will be available in leafviz through the Summary tab). [%default]"))
)
parsed_args <- parse_args(option_parser, positional_arguments = 4)
@@ -34,6 +27,12 @@ results_file = parsed_args$options$output
groups_file <- parsed_args$options$meta_data_file
FDR_limit <- parsed_args$options$FDR
+require(leafcutter)
+library(data.table)
+library(stringr)
+library(dplyr)
+library(magrittr)
+
cat("Preparing for visualisation\n")
# annotation