This repository contains the scripts used to process 16S and ITS data from soils collected at the Savannah River Site (SRS) in May 2021. Scripts were written by Claire Winfrey, with substantial input on analyses from Julian Resasco and Noah Fierer. Below is a description of each script, listed in the order that scripts should be run to replicate results. Analyses on 16S and ITS data were done in separate scripts, in parallel.
16S
- SRS_16S_Bioinformatics.R- This script demultiplexes raw 16S rDNA reads from the MiSeq, uses dada2 to trim, denoise, dereplicate, and merge reads, and assigns bacterial/archaeal taxonomy (using SILVA version 132).
- 16SExploratoryDataAnalysisAug2021.R - This script takes a first look at 16S MiSeq data. I look at sequence counts, rarefy, examine top classes and phyla, and look at ordinations bsaed on different EUs, different plant communities, etc. I repeat the ordinations and examination of top classes and phyla before and after removing rare species. I also investigate samples that seem to be outliers in ordination space. Finally, this script creates the cleaned up dataset used in all downstream scripts and analyses.
- UbiquityMedianSetup.R- This script applies a ubiquity filter to the "trimmed" dataset made in 16SExploratoryDataAnalysisAug2021.R.
- DissimilaritiesAcrosTransectsMedianAbund.R -- This script calculates and plots Bray-Curtis dissimilarities across the transects using the median ASV abundance object (medianEU.ps) made in UbiquityMedianSetup.R. A dbRDA and model selection are performed here, but are probably not the best to use since they do not incorporate the transect-level difference (i.e. are just based on the median).
- DiffAbundProkaryotes.R - This script performs a differential abundance analysis using DESeq2 package to find patch and forest specialist ASVs for ONLY bacterial and archaeal data (does so for dataset with all EUs and for all EUs EXCEPT 53N). Importantly, this script uses the post ubiquity dataset created in UbiquityMedianSetup.R (with and without 53N), but with 1 count added to each of the ASV abundance counts (necessary for the ITS DESeq analysis). This script then makes a stacked barplot showing the number of ASVs in each phyla of the postUbiquity dataset, and whether it is an indicator species for the patch, the matrix, or neither. Finally, this script makes plots comparing the ubiquity of different specialist groups.
- ProkaryotesAbundanceZScores.R - Using the differentially abundant taxa identified in DiffAbundProkaryotes.R, this script transforms the abundance of each of these ASVs into a z-score (where the mean and standard deviation for calculation of that ASV's z-score are drawn from within a given EU but at each point along the transect). Then, using these z-scores, the script finds the logistic fit and associated parameters of the curve of each ASV as it goes over the transect. Finally, this script makes visualizations of these results.
- PostUbiqGraphics_16S.R- makes NMDS ordinations and Bray-Curtis dissimilarities along the transect for dataset created in UbiquityMedianSetup.R (i.e. after removing rare ASVs and not very ubiquitous ones). In addition, creates stacked barcharts of the relative abudnance of top phyla and families. Finally, performs a PERMANOVA (adonis function) and then a dbRDA (ultimately in manuscript) to test if samples segregate based on habitat type.
ITS
- SRS_ITS_Bioinformatics.R- This script demultiplexes raw ITS rDNA reads from the MiSeq, uses dada2 to trim, denoise, dereplicate, and merge reads, and assigns bacterial/archaeal taxonomy (using UNITE database, Oct. 5, 2021 release).
- ITS_ExploratoryDataAnalysis.R - This script explores the ITS sequencing data, roughly doing the same thing for ITS data as 16SExploratoryDataAnalysisAug2021.R. However, unlike in the 16S script, no differential abundance analysis is done at this stage to explore patch versus matrix.
- ITS_UbiquityMedianSetup.R- this script applies a ubiquity filter to the "trimmed" dataset made in ITS_ExploratoryDataAnalysis.R. Also creates a dataset where ASV abundance is now median (values are across 4 transects at a given meter location within each EU). Finally, this script formats an ASV table for compatibility with FUNGuild.
- fungiDissimilaritiesAcrosTransectsMedianAbund.R -- This script calculates and plots Bray-Curtis dissimilarities across the transects for fungi using the median ASV abundance object (ITS_medianEU.ps) made in ITS_UbiquityMedianSetup.R. The script is almost identifical to DissimilaritiesAcrosTransectsMedianAbund.R. A dbRDA and model selection are performed here, but are probably not the best to use since they do not incorporate the transect-level difference (i.e. are just based on the median
- DiffAbundFungi.R-- this script performs a differential abundance test to see which ASVs differ between the forested matrix and the savanna patch, and saves these results (which are loaded in the following script). It also makes visualizations of these results. Basically identical to DiffAbundProkaryotes.R except with fungi.
- FungiAbundanceZscores.R - same as ProkaryotesAbundanceZScores.R except with fungi.
- PostUbiqGraphic_ITS.R - makes NMDS ordinations and Bray-Curtis dissimilarities along the transect for dataset created in ITS_UbiquityMedianSetup.R (i.e. after removing rare ASVs and not very ubiquitous ones). In addition, creates stacked barcharts of the relative abudnance of top phyla and families. Finally, performs a PERMANOVA (adonis function) and then a dbRDA (ultimately in manuscript) to test if samples segregate based on habitat type.
- FUNGuildExploraton.R - This script prepares ASV dataset (ASV table made in ITS_UbiquityMedianSetup.R, where ASVs hve been rarefied and ubiquity and rarity filters applied) for FUNGuild and then explores results of FUNGuild and makes stacked barplot visualizations. Run up to line 72, then run FUNGuildStepByStep, then return to this script to explore results.
- FUNGuildStepByStep text file describes the specific steps taken to run FUNGuild, as well as provides Python code run in command line. Should be run after running through line 72 in FUNGuildExploraton.R.
16S and ITS together or other miscellaneous
- EnvDataPlots.R -- this script makes plots of mean pH, canopy cover, and vegetation for all of the points along the transects sampled in the dataset. In addition, it looks at soil analysis for transect end points.
- DissimAcrossTransectDiffAbundOnly.R -- this script looks at Bray-Curtis dissimilarities across the transect within each EU, focusing only on ASVs found to be differentially abundant (in DiffAbundProkaryotes.R and DiffAbundFungi.R).
- AbundanceZscores_no53N.R -- This script calculates Z-scores, with EU 53N removed from analysis.
- AbundanceDropOffAUCsAll.R -- This script gets the area under the curve for each indicator ASV in each EU, and in each one of these ASV/EU combinations, finds the meter where the AUC is half that of the maximum.
Scripts with additional analyses that probably won't make it into paper
- EdgeEffectsAllSitesFUNGI.R-- see description above for EdgeEdgebyASV_allSites.R, but this script uses the ITS data. Was replaced with DiffAbundFungi.R.
- EdgeEdgebyASV_allSites.R-- was replaced with DiffAbundProkaryotes.R in work flow of data analysis.