Last updated: 2024-01-24
Checks: 7 0
Knit directory:
locust-comparative-genomics/
This reproducible R Markdown analysis was created with workflowr (version 1.7.1). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
The command set.seed(20221025) was run prior to running
the code in the R Markdown file. Setting a seed ensures that any results
that rely on randomness, e.g. subsampling or permutations, are
reproducible.
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.
Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.
Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.
The results in this page were generated with repository version 838b0a1. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.
Note that you need to be careful to ensure that all relevant files for
the analysis have been committed to Git prior to generating the results
(you can use wflow_publish or
wflow_git_commit). workflowr only checks the R Markdown
file, but you know if there are other scripts or data files that it
depends on. Below is the status of the Git repository when the results
were generated:
Ignored files:
Ignored: .DS_Store
Ignored: .RData
Ignored: .Rhistory
Ignored: .Rproj.user/
Ignored: analysis/.DS_Store
Ignored: analysis/.Rhistory
Ignored: analysis/site_libs/
Ignored: data/.DS_Store
Ignored: data/.Rhistory
Ignored: figures/
Unstaged changes:
Modified: analysis/1_team-science.Rmd
Modified: analysis/1_workenvs-prep.Rmd
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
These are the previous versions of the repository in which changes were
made to the R Markdown (analysis/1_seq-data.Rmd) and HTML
(docs/1_seq-data.html) files. If you’ve configured a remote
Git repository (see ?wflow_git_remote), click on the
hyperlinks in the table below to view the files as they were in that
past version.
| File | Version | Author | Date | Message |
|---|---|---|---|---|
| Rmd | 838b0a1 | Maeva A. TECHER | 2024-01-24 | wflow_publish("analysis/1_seq-data.Rmd") |
| Rmd | b195f57 | Maeva A. TECHER | 2024-01-23 | refresh |
| html | b195f57 | Maeva A. TECHER | 2024-01-23 | refresh |
| html | 612efb5 | Maeva A. TECHER | 2023-12-18 | Build site. |
| Rmd | 53877fa | Maeva A. TECHER | 2023-12-18 | add pages |
Load R libraries (install first from CRAN or Bioconductor)
library("knitr")
library("rmdformats")
library("tidyverse")
library("DT") # for making interactive search table
library("plotly") # for interactive plots
library("ggthemes") # for theme_calc
library("reshape2")
## Global options
options(max.print="10000")
knitr::opts_chunk$set(
echo = TRUE,
message = FALSE,
warning = FALSE,
cache = FALSE,
comment = FALSE,
prompt = FALSE,
tidy = TRUE
)
opts_knit$set(width=75)The Behavioral Plasticity Research Institute (BPRI) successfully generated high-quality, chromosome-length genomes for six species of grasshoppers and locusts within the Schistocerca genus. To achieve this, the genome assembly team utilized a hybrid approach, combining long-read sequencing with HiFi PacBio technology and short-read sequencing with Hi-C. Following the finalization of the genome assemblies, the annotation process was conducted by the Eukaryotic Annotation Pipeline for RefSeq by NCBI.
We will use the assemblies that have received accession numbers in RefSeq and are associated with .gtf and .gff files.
here is an an example of how to download the genome assembly for Schistocerca piceifrons from NCBI RefSeq:
## Downloading the genome sequence, primary assembly fasta file with RefSeq contig accessions
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/021/461/385/GCF_021461385.2_iqSchPice1.1/GCF_021461385.2_iqSchPice1.1_genomic.fna.gz
## Downloading the comprehensive annotation file in .gtf format
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/021/461/385/GCF_021461385.2_iqSchPice1.1/GCF_021461385.2_iqSchPice1.1_genomic.gtf.gzThe reads that can be mapped on these sequences can be either from public sequences with accession on NCBI, DDBJ or generated de novo and yet to be released.
We conducted a search in the NCBI database for RNA SRA associated with the Acrididae, Schistocerca genus and found 160 accessions were available. All RNA-seq data for species including S. americana, S. piceifrons, S. serialis cubense, S. cancellata, S. nitens were originally generated in the Song Laboratory of Insect Systematics and Evolution at Texas A&M University and under the same conditions (Foquet et al. 2021).
Using the Run Selector from NCBI, we can easily download a metadata table that allows us to visualize how the accessions are distributed per species. We will use this metadata table for analysis of the RNAseq data from the previously published dataset.
We compile a list of accessions from the “Run” category for
each species and then use SRA-toolkit from NCBI. First, we
create an empty directory named ncbi to download each SRA.
This is where SRA Toolkit will dump the prefetched SRA
files in the .sra format.
ml purge
ml GCC/10.2.0 OpenMPI/4.0.5 SRA-Toolkit/2.10.9
vdb-config --interactiveOnce in the vdb-config interactive mode, select cache,
choose, then use [ .. ], to enter
/home/USERNAME/PATH/ncbi one directory at a time
prefetch --option-file SraAccList.txt
cat SraAccList.txt | xargs fasterq-dump --split-3 --outdir "/your-directory/for-fastq"Clean-up the ncbi directory and move the fastq.gz file
(rename if wanted).
another option
for x in *.sra ; do fasterq-dump --split-files $x ; mv *.fastq ../../paired_end_piceifrons/; doneWe sequenced new whole transcriptomes from the desert locust S. gregaria under two different density conditions. While these new transcriptomes originated from specimens reared in similar conditions, we made slight changes in the RNA extraction, library preparation and sequencing processes. Specifically, we used the Illumina Stranded Total RNA with RiboZero depletion kit and sequenced on a NovaSeq SP flow cell at TxGen, aiming for a targeted yield of 40 millions reads per library.
For the sequencing of many transcriptomes at the same time, it is common practice to run pooled libraries across multiples lanes, depending on the required read count. To handle this, I used the following loop to merge several files from different lanes:
#!/bin/bash
# add the print $3 for TxGen reads as it is for samples S
for i in $(ls -1 *R1*.gz | awk -F '_' '{print $1"_"$2"_"$3}' | sort | uniq)
do echo $i
echo "Merging R1 ${i}"
cat "$i"_L00*_R1_001.fastq.gz > "$i"_MERGE_R1_001.fastq.gz
echo "Merging R2 ${i}"
cat "$i"_L00*_R2_001.fastq.gz > "$i"_MERGE_R2_001.fastq.gz
done;We need to create a .csv file containing as much
information as possible for each sample/file name (e.g., Sample_ID,
Species, Sex, RearingCondition). An interactive and searchable table
with this type of information is available below and can even be
downloaded directly.
NB: Throughout our analysis, we will complete this metadata file by adding other stats related to sequencing and mapping.
# Load our SRA metadata table
metaseq <- read_table("data/metadata/RNAseq_modified_METADATA2022.txt", col_names = TRUE,
guess_max = 5000)
## Create an interactive search table
metaseq %>%
datatable(extensions = "Buttons", options = list(dom = "Blfrtip", buttons = c("copy",
"csv", "excel"), lengthMenu = list(c(10, 20, 50, 100, 200, -1), c(10, 20,
50, 100, 200, "All"))))
sessionInfo()FALSE R version 4.3.1 (2023-06-16)
FALSE Platform: x86_64-apple-darwin20 (64-bit)
FALSE Running under: macOS Sonoma 14.1.2
FALSE
FALSE Matrix products: default
FALSE BLAS: /Library/Frameworks/R.framework/Versions/4.3-x86_64/Resources/lib/libRblas.0.dylib
FALSE LAPACK: /Library/Frameworks/R.framework/Versions/4.3-x86_64/Resources/lib/libRlapack.dylib; LAPACK version 3.11.0
FALSE
FALSE locale:
FALSE [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
FALSE
FALSE time zone: America/Chicago
FALSE tzcode source: internal
FALSE
FALSE attached base packages:
FALSE [1] stats graphics grDevices utils datasets methods base
FALSE
FALSE other attached packages:
FALSE [1] reshape2_1.4.4 ggthemes_5.0.0 plotly_4.10.4 DT_0.31
FALSE [5] lubridate_1.9.3 forcats_1.0.0 stringr_1.5.1 dplyr_1.1.4
FALSE [9] purrr_1.0.2 readr_2.1.5 tidyr_1.3.0 tibble_3.2.1
FALSE [13] ggplot2_3.4.4 tidyverse_2.0.0 rmdformats_1.0.4 knitr_1.45
FALSE [17] workflowr_1.7.1
FALSE
FALSE loaded via a namespace (and not attached):
FALSE [1] gtable_0.3.4 xfun_0.41 bslib_0.6.1 htmlwidgets_1.6.4
FALSE [5] processx_3.8.3 callr_3.7.3 tzdb_0.4.0 crosstalk_1.2.1
FALSE [9] vctrs_0.6.5 tools_4.3.1 ps_1.7.5 generics_0.1.3
FALSE [13] fansi_1.0.6 pkgconfig_2.0.3 data.table_1.14.10 lifecycle_1.0.4
FALSE [17] compiler_4.3.1 git2r_0.33.0 munsell_0.5.0 getPass_0.2-4
FALSE [21] httpuv_1.6.13 htmltools_0.5.7 sass_0.4.8 yaml_2.3.8
FALSE [25] lazyeval_0.2.2 crayon_1.5.2 later_1.3.2 pillar_1.9.0
FALSE [29] jquerylib_0.1.4 whisker_0.4.1 ellipsis_0.3.2 cachem_1.0.8
FALSE [33] tidyselect_1.2.0 digest_0.6.34 stringi_1.8.3 bookdown_0.37
FALSE [37] rprojroot_2.0.4 fastmap_1.1.1 grid_4.3.1 colorspace_2.1-0
FALSE [41] cli_3.6.2 magrittr_2.0.3 utf8_1.2.4 withr_2.5.2
FALSE [45] scales_1.3.0 promises_1.2.1 timechange_0.2.0 rmarkdown_2.25
FALSE [49] httr_1.4.7 hms_1.1.3 evaluate_0.23 viridisLite_0.4.2
FALSE [53] rlang_1.1.3 Rcpp_1.0.12 glue_1.7.0 formatR_1.14
FALSE [57] rstudioapi_0.15.0 jsonlite_1.8.8 R6_2.5.1 plyr_1.8.9
FALSE [61] fs_1.6.3