This dataset comes from the following paper by researchers at JCU and the Australian Institute of Marine Science (AIMS);
Destabilization of mutualistic interactions shapes the early heat stress response of the coral holobiont
The full paper can be found at this link: Marangon et al 2025
When corals are exposed to prolonged heat stress they bleach, meaning they lose their symbiotic microalgal partners. These symbionts (unicellular algae) give corals most of their color, so when they are lost the corals turn white. However, this description only considers two components of the coral holobiont. A third component that is often overlooked are the bacteria and archaea associated with the coral. This experiment addresses that gap by examining how bacterial and archaeal (via 16 microbial profiling) communities in a coral change with heat stress.
Fragments of the coral Porites lutea were cut from three parent colonies (A, B and C) on the Great Barrier Reef. Fragments of the same colony are genetically identical, so colony is effectively coral genotype. The fragments were spread across six aquarium tanks so that every tank held fragments of all three colonies. Three tanks were kept at ambient temperature (about 28.5 °C). The other three were heated slowly, from 28.2 °C to 32 °C over 56 days, until the corals had experienced 8 degree heating weeks (DHW). This is the level of heat stress at which severe bleaching is expected on reefs.
At each sampling time one fragment of each colony was removed from each tank, and seawater from each tank was filtered to sample the bacteria living in the water around the corals. The V4 region of the bacterial 16S rRNA gene was sequenced from each sample. This dataset includes samples taken at the start of heating (T0, day 0) and at the end (T5, day 56, 8 DHW in the heated tanks).
| Sample type | T0 | T5 |
|---|---|---|
| Corals, ambient tanks | 9 | 9 |
| Corals, heated tanks | 9 | 9 |
| Seawater, ambient tanks | 3 | 3 |
| Seawater, heated tanks | 3 | 3 |
Four additional samples are present. These are baseline samples, taken before the experiment was started (before fragmentation, after acclimation). There is one sample each for colony B and C and two samples from colony A.
This assignment starts from raw reads. You won’t need to demultiplex but you will need to perform all the steps after demultiplexing in the 16S microbiome workflow. This includes long-running steps such as removing primers and denoising the reads yourself. Start early because there will be times when you need to wait for jobs to finish on the server.
First create a directory raw_data for storing the original data for this project and another directory called cache for storing your own analysis outputs
mkdir raw_data
mkdir cache
The data for this assignment is located at the path /pvol/data/porites on the class rstudio server. Run the following commands to create a symbolic link to this in the raw_data directory of your assignment
cd raw_data
ln -s /pvol/data/porites .
The following files are present within the porites data directory;
| File Name | Purpose |
|---|---|
fastq/ |
Raw reads. There are two files per sample, <sample-id>_R1.fastq.gz and <sample-id>_R2.fastq.gz. These are Illumina MiSeq 2 x 250 bp paired-end reads in compressed fastq format. |
manifest.tsv |
A manifest listing the fastq files for each sample, for importing the reads with qiime tools import |
sample-metadata.tsv |
Sample Metadata |
classifier/silva-138-99-515-806-nb-classifier-slim.qza |
Naive Bayesian classifier trained using SILVA 138. It was trained on the region amplified by the 515F and 806R primers used in this study, and it has been modified to use less memory. It is appropriate for use with the --i-classifier argument required by qiime feature-classifier classify-sklearn. |
The most important columns in the sample metadata file are;
AMB or HEAT), the time point, the colony and the tank. For example HEAT-T5-C-t24 is the colony C fragment taken from heated tank 24 at T5. Seawater samples start with SW and baseline corals with BASE.heat_T5 (seawater groups start with seawater_)t2, t9, t23 (ambient) or t7, t13, t24 (heat); empty for baseline coralsTemperature was applied to whole tanks, so each treatment has only three true replicates (its tanks), not nine. Think about how tank and colony should be accounted for in your statistical tests (for example in the formula given to qiime diversity adonis).
Unlike the metagenomics walkthrough exercise, this time we have pre-demultiplexed data (much more common these days) so the workflow for importing data into a qiime object is slightly different. In this case we will supply a file manifest.tsv with all the file paths for the raw fastq files. QIIME will read all these files and import them to a single object. This is the command you’ll need
qiime tools import \
--type 'SampleData[PairedEndSequencesWithQuality]' \
--input-path raw_data/porites/manifest.tsv \
--input-format PairedEndFastqManifestPhred33V2 \
--output-path cache/demux.qza
This data will need trimming to remove adapters. You should use the cutadapt plugin in qiime for this as follows;
qiime cutadapt trim-paired \
--i-demultiplexed-sequences cache/demux.qza \
--p-front-f GTGYCAGCMGCCGCGGTAA \
--p-front-r GGACTACNVGGGTWTCTAAT \
--p-discard-untrimmed \
--p-cores 4 \
--o-trimmed-sequences cache/trimmed.qza
After this finishes I recommend you create a visualisation for the trimmed reads as follows
qiime demux summarize \
--i-data cache/trimmed.qza \
--o-visualization cache/trimmed.qzv
Use trimmed.qzv to look at read quality before you choose truncation lengths for denoising.
Long-running commands should be run as slurm jobs, as you did in the metagenomics walkthrough. For these you’ll need to make sure you define the qiime alias as the top of your script. Here is an example script which performs the import command shown above.
#!/bin/sh
#SBATCH --time=60
#SBATCH --ntasks=1 --cpus-per-task=4 --mem=2gb
alias qiime='apptainer run -B /pvol/:/pvol --env "MPLCONFIGDIR=$HOME/matplotlib" /pvol/data/sif/qiime.sif qiime'
qiime tools import \
--type 'SampleData[PairedEndSequencesWithQuality]' \
--input-path raw_data/porites/manifest.tsv \
--input-format PairedEndFastqManifestPhred33V2 \
--output-path cache/demux.qza
Although some qiime commands should run quickly there are four that take a long time and require a lot of memory. This table shows the approximate run times and memory. When you perform each of these steps you MUST do it using a slurm script.
These are approximate times and peak memory for the slow steps with 4 cores. They will be longer when the server is busy, so ask slurm for more time than this.
| Step | Time | Memory |
|---|---|---|
| Import, primer removal and summary | 15 min | 1 GB |
qiime dada2 denoise-paired |
1.25 h | 5 GB |
qiime feature-classifier classify-sklearn |
2 h | 10 GB |
qiime phylogeny align-to-tree-mafft-fasttree |
15 min | 1 GB |
percentage of input merged in the denoising stats to see if you have lost too many reads (reads that don’t overlap will be discarded).--p-reads-per-batch 250 with classify-sklearn. Without it, every ASV goes into a single batch and classification needs much more memory. There is a slurm script for classification below.qiime taxa filter-table before calculating diversity.This slurm script classifies your ASVs with 4 cores in about 2 hours;
#!/bin/sh
#SBATCH --time=4:00:00
#SBATCH --ntasks=1 --cpus-per-task=4 --mem=12gb
alias qiime='apptainer run -B /pvol/:/pvol --env "MPLCONFIGDIR=$HOME/matplotlib" /pvol/data/sif/qiime.sif qiime'
qiime feature-classifier classify-sklearn \
--i-classifier raw_data/porites/classifier/silva-138-99-515-806-nb-classifier-slim.qza \
--i-reads cache/rep-seqs.qza \
--p-n-jobs 4 \
--p-reads-per-batch 250 \
--o-classification cache/taxonomy.qza
As this is a 16S Metabarcoding dataset you should complete the following basic analyses in qiime2
qiime dada2 denoise-pairedqiime feature-table summarize to guide a decision on the rarefaction thresholdVisualisation and statistical tests should focus on the question(s);