Experimental design

To simulate the effects of overfishing, nutrient loading or the combination of these stressors, we conducted a 3-year field experiment. Four pairs of 9-m2 plots were established. One member of each of these pairs was enriched with nitrogen and phosphorous, while the other remained at ambient nutrient levels (Supplementary Fig. 1). These plots were >10 m from each other in all cases. Each 9-m2 plot was delineated into nine 1-m2 subplots with metal nails driven into the reef at the corners and centre of each plot. The locations of the plots were selected such that initial variation in rugosity and algal cover within each subplot was minimal. Within each plot, two randomly selected subplots were enclosed with herbivore exclosures, while two other random subplots were selected as exclosure controls. Exclosure controls were fitted with open-topped exclosures. These controls allowed access by herbivorous fishes, but acted as controls for other potential artifacts of the cages.

All exclosures were made of plastic-coated wire mesh with 2.5-cm diameter holes. This diameter mesh generally excludes most fishes >10 cm total length. Smaller or juvenile herbivorous fishes are able to enter the exclosures, but these smaller herbivores generally contribute little to overall grazing rates on reefs and have minimal impacts on the algal communities59. In addition, access by smaller herbivores reflects patterns seen under intensive fishing, in which larger fish species are preferentially harvested while leaving smaller size classes of fish60,61. We scrubbed both exclosures and exclosure controls every 4–6 weeks to remove fouling organisms.

Nutrient pollution was simulated using slow-release fertilizer diffusers applied to each nutrient enrichment plot. Each diffuser was a 15-cm diameter PVC tube, perforated with six 1.5 cm holes. The open ends of the PVC tube were wrapped in fine plastic mesh to keep fertilizer pellets inside, but allow diffusion of soluble nutrients. A total of 175 g of Osmocote (19–6–12, N–P–K) slow-release fertilizer was loaded into each diffuser. PVC enrichment tubes were attached to each metal nail within the 9-m2 enrichment plots for a total of 25 enrichment tubes per enrichment plot. Nutrients were replaced every 30–40 days to ensure continued delivery of N and P. Previous studies have shown Osmocote delivery using similar methods to be an effective way of enriching water column nutrients in benthic systems (for example, ref. 62).

Confirmation of the efficacy of nutrient enrichment

Nitrogen and phosphorus levels were assessed in the water column above each enrichment and control plot as in ref. 52. In July of 2009, 2010 and 2011, divers used 60-ml syringes to slowly draw water from ∼3 cm above the benthos in either control or enriched plots. Samples were taken ∼30 days after nutrient diffuser deployment to ensure that enrichment occurred across the duration of diffuser deployment. Immediately after collection, samples were filtered (GF/F) into acid-washed bottles, placed on ice, returned to the laboratory and frozen until analysed. Dissolved inorganic nitrogen (DIN=ammonium and nitrite+nitrate) and soluble reactive phosphorus concentrations were determined via autoanalyser. We also assessed nutrient enrichment efficiency by analysing tissue carbon:nitrogen (C:N) levels in the common alga Dictyota menstrualis. The nutrient content of macroalgae such as D. menstrualis reflects ambient nutrient conditions over longer time scales (that is, weeks to months) than ambient water column nutrients (ref. 63 and references therein). We collected D. menstrualis from both enriched and control treatments during the same months as water samples. Tissues were dried at 60 °C, ground to a powder and analysed for %C and %N content with a CHN Carlo-Erba elemental analyzer (NA1500). Nutrient data from both water and algal tissue for each replicate were averaged across summers for statistical analysis via analysis of variance (ANOVA).

Quantification of the herbivorous fish community

Periodically throughout the study, we used 30 × 2 m belt transects (n=8) to quantify the density of different herbivorous fishes. Divers slowly swam the length of each transect counting individuals of the different herbivorous fishes in the genera Sparisoma, Scarus, Acanthurus and Kyphosus. Fishes were identified to species and their length was estimated to the nearest centimetre. We used published length:weight relationships64 to convert fish densities into herbivore biomass. We analysed these data with one-factor ANOVA examining potential differences in herbivore biomass over time.

Quantification of benthic cover

At least once every season (for example, spring, summer, fall, winter at 12–14 week intervals), we visually quantified benthic cover within four, 50 × 50 cm quadrats in each of the 1 m2 treatment areas. These quadrats were divided into 49 points, and benthic organisms under each point were identified to species or genus. Algae that are challenging to identify taxonomically under field conditions (for example, crustose coralline algae and filamentous algae) were classified into algal functional groups. Filamentous algae were classified into short algal turf (<0.5 cm in height) or algal turf (>0.5 cm in height) given that taller, thicker algal turf can often be deleterious to coral health and growth10.

Benthic cover was quantified in June 2009 1 week before treatments were initiated to provide a baseline from which to assess changes in algal abundance and community structure. No significant differences among treatments in algal abundance could be detected at the beginning of the experiment (see initial time points in Fig. 1a,b), as expected given random assignment of subplots to treatment conditions. Further, during the summer of each year (2009–2012) when algal cover was often at its highest, we also surveyed open areas of reef (areas that did not have three-sided exclosure controls) within the 9-m2 plots to assess whether the exclosure controls had any effect on algal abundance or community composition. We did not detect any differences in algal abundance or community composition between the open unmanipulated areas and exclosure controls (Supplementary Data 1).

Coral tissue growth or loss analyses

At the beginning of the experiment, we mapped each coral colony in the experimental plots that were >2 cm in diameter and took close-up photographs of these corals in situ. Subsequently, we photographed each of these corals every ∼16 weeks throughout the experiment for a photographic record of changes in coral colony health. In each picture a ruler or object of known size was placed next to the coral to provide scale. In total, we tracked the fate of 226 individual corals spread across each of the treatments for over 3 years. The most common corals were Porites porites (41.1 % of corals), Agaricia spp. (17.7 % of corals), Siderastrea siderea (15.5 % of corals) and P. astreoides (11.5 % of corals).

These corals allowed us to evaluate the impact of the different treatments on coral growth or tissue loss across the time course of the experiment. We scored growth or tissue loss on a 12-point scale, with bins corresponding to amounts of tissue loss that could be readily observed in photographs (for example, −2=10–25% tissue loss). We scored the tissue loss or gain of each coral over the course of the experiment on the following scale: −6=100% tissue loss, −5=75–90% loss, −4=50–75% loss, −3=25–50% loss, −2=10–25% loss, −1=0–10% loss, 0=0% loss/gain, 1=0–10% gain, 2=10–25% gain, 3=25–50% gain, 4=50–75% gain, 5=75–90% gain and 6>100% gain. We then converted these scores to mean loss/gain by averaging the range corresponding to that score. For example, a coral with a −3 score would be converted to a −37% tissue loss value. Only nine corals grew >100% (score=6) over the course of the experiment. For these corals, we estimated the growth for each coral at 100–500% at 50% intervals (for example, 100, 150, 200% and so on). Statistical analyses were conducted based on the raw tissue gain/loss scores, but converted to percentages in the presentation for ease of interpretation. Further, at each time point we scored each coral for: (1) algal competition as measured by direct contact with algal competitors (and the identification of that algal competitor), (2) the presence of overlying sediment on the coral, (3) predation scars from parrotfishes and invertebrate corallivores (only the former were observed at appreciable levels), and (4) signs of bleaching or disease. The primary coral disease observed was DSS (see ref. 52 for additional discussion).

Statistical analysis of benthic cover and coral health

Algal cover was analysed using mixed effects models to determine whether the response variable differed among enrichment treatment, herbivore treatment and seasons, as well as whether there was an enrichment × herbivore × season interaction (Supplementary Data 1 sheet c). The nested, split-plot design of the experiment was incorporated into the model by nesting replicates of the exclosures and exclosure controls within ambient or nutrient-enriched plots. We analysed cover for important species or functional groups, as well as for overall upright algal cover, which is a proxy for the competitive environment of corals. Upright algal cover included all macroalgae and tall filamentous turf, but excluded crustose coralline algae and short filamentous turf as these two functional groups are relatively benign for corals10. In addition, we assessed how herbivore exclusion, nutrient pollution and season impacted algal community structure via PCoA of Bray–Curtis divergences, as well as permutational MANOVA (PERMANOVA).

Per cent coral mortality per treatment and coral tissue loss were analysed using similar mixed models to algal cover. For growth measures, corals were nested within ambient or enriched plots, but we did not incorporate season as we only analysed change in tissue for corals at the end of the experiment. We calculated tissue loss statistics either excluding or including corals that suffered total colony mortality. Corals that died suffered total colony morality and therefore 100% loss of live tissue area. Including these corals in coral growth analyses resulted in non-normal distributions that could not be corrected via transformations. Therefore, we analysed coral growth both excluding the corals that died, which satisfied normality requirements for the analyses, and including the corals that died. Both analyses produced relatively consistent results (Supplementary Data 5), with the exception that only the interaction of herbivory × nutrient loading was significant in Porites corals (rather than each factor also being individually significant) when total colony mortality was excluded. We used a χ2-test to determine if coral mortality was higher or lower than expected across different seasons given the null hypothesis that coral mortality would be distributed evenly across seasons (25% of total mortality per season).

To assess the impact of algal competition and parrotfish predation on coral mortality for Porites spp., (the corals with the highest mortality rates), we used Fisher’s exact test to test for differences in the proportion alive versus dead corals across enrichment and/or herbivore exclusion treatments. We also used Fisher’s exact tests to determine how factors such as competition with different types of algae (that is, corals with or without algal contact) led to loss/gain of coral tissue or coral death/survival.

Mixed effects models were run using the ‘nlme’ package and post hoc comparisons conducted using the ‘multcomp’ package in R v3.0.0. Fisher’s exact tests on contingency tables were run using JMP software (SAS). PCoA of Bray–Curtis divergences and PERMANOVA of benthic community data were conducted in QIIME 1.8.

Field sampling of coral mucus for microbial analysis

Coral mucus microbial communities were studied in depth for three coral genera common to the study region and most abundant within the plots: Siderastrea (Siderastrea siderea only), Porites and Agaricia. We used 16S rRNA gene surveys to study the microbial changes in the coral mucus across the course of treatment, focusing especially on whether detectable microbial changes accompanied specific routes to mortality or tissue loss in corals. (For metadata and sequencing statistics, see Supplementary Data 2; for OTU table in BIOM format, see Supplementary Data 7). We focused on mucus communities because these are thought to provide a barrier against invasion by opportunistic pathogens, and can be sampled non-destructively from the same individuals over time. We deemed other methods, such as collecting live coral tissue, too harmful and invasive to the coral for our goal of monitoring the coral microbiome over the long term.

Coral-associated bacteria and archaea were collected using sterile syringe removal of the coral surface mucus layer on SCUBA. A sterile syringe was used to first agitate the coral and then remove 10 ml of mucus using negative pressure. Once on the surface, the sample was placed in a sterile 15 ml conical tube, frozen on dry ice for the return trip and kept frozen until processing. At a number of time points, 15 ml water sample controls were collected from >1 m above the reef and treated identically to mucus samples in downstream processing. In this study, we did not attempt to evaluate changes in Symbiodinium abundance or taxonomy, as it has not yet been established that mucosal abundances of Symbiodinium types reflect abundances within tissues, and destructive tissue sampling would preclude time-series analysis of the same coral colonies.

Quantification of sub-bleaching thermal stress

We calculated cumulative thermal stress using the DHW metric of National oceanic and atmospheric administration (NOAA’s) Coral Reef Watch47. These calculations rely on a climatological baseline for the mean temperature of the warmest month known as the MMM. We calculated the MMM for our site using the NOAA Pathfinder v5.2 data set, the US official climate data record for sea-surface temperature45, and found that the MMM for our site was 29.26 °C (Supplementary Data 6). Data used spanned 1982–2008, and excluded study dates. Temperatures that exceed the MMM by +1 °C are generally regarded as constituting thermal stress47. We calculated an MMM +1 °C threshold of 30.26 °C for our site. Temperatures above the 30.26 °C threshold for accumulation of DHWs at our site only occurred during the warm seasons (3 months centred on the warmest month, August).

DHWs represent the extent to which temperatures exceeded the MMM +1 °C in a given season. Temperatures were above 30.26 °C for 7 weeks (total) during all sampling years, resulting in the accumulation of <1 DHW during each of the study years (Supplementary Fig. 10a). In predictions of coral bleaching, accumulation of four DHWs is often associated with minor to moderate bleaching5. We saw little bleaching within experimental plots, consistent with sub-bleaching levels of thermal stress. Sub-bleaching thermal stress may nonetheless negatively affect reef health if it increases the abundance and/or virulence of bacterial pathogens. For example, one recent model15 proposes that bacterial pathogens may become problematic for corals at temperatures, exceeding the MMM (a lower threshold than the MMM + 1 °C threshold used in coral bleaching studies). We therefore related our microbiological data to both the MMM +1 °C for coral thermal stress during coral bleaching, and the MMM threshold proposed for bacterial pathogens (Supplementary Fig. 10b).

Coral mucus microbiome sample processing and sequencing

In the laboratory, coral mucus samples were thawed, centrifuged and supernatant decanted. DNA was purified using an organic extraction as previously described26. After DNA extraction, microbial 16S amplicon libraries were generated using the primers 515F and 806R, both with added 454 sequencing adaptors and with Golay barcodes added to the reverse primer. Triplicate 25 μl reactions were conducted using the GoTaq Flexi system from Promega (Madison, WI, USA) with the following conditions per reaction: 1 × clear buffer, 1 mM dNTPs, 5 mM MgCl, 1 μM of each primer, 1u Taq polymerase and 1 μl of extracted DNA template. Thermocycling was conducted as follows: 1 cycle of 94 °C for 3 min; 35 cycles of 94 °C for 45 s, 50 °C for 60 s and 72 °C for 90 s; and 1 cycle of 72 °C for 10 min. Amplification success was checked on a 1.5% agarose gel, and successful triplicate reactions were pooled and cleaned using AMPure magnetic beads from Agencourt. Before sequencing, libraries were quantified using a Qubit dsDNA HS kit from Invitrogen and then pooled into equimolar ratios. The pooled library was checked for amplicon length and purity on an Agilent Bioanalyzer 2100 and then sequenced on a 454 Roche pyrosequencer (GSJunior platform) at the Oregon State University’s Center for Genome Research and Biocomputing Core Laboratories.

Microbial community data quality control

The QIIME (v.1.8) software pipeline65 was used for quality control, selection of operational taxonomic units and analyses of community diversity (Supplementary Data 2). Sequence libraries were demultiplexed, and sequences with quality scores less than a mean of 35 were removed. Error-correcting barcodes were used to detect and recover sequences whose barcode sequence had exactly one sequencing error. Barcode sequences with two or more errors were removed. Sequences were clustered into operational taxonomic units (OTUs), at a 97% 16S rRNA gene identity threshold using USEARCH 6.1.544, and the subsampled open-reference OTU-picking protocol in QIIME v.1.8 (ref. 65, using greengenes 13_8 as the reference66. This OTU-picking protocol clusters all reads, but assigns reference ids to OTUs in greengenes, which can be useful in comparisons across studies. Chimeric sequences were removed with UCHIME. OTUs represented in the overall analysis by only a single count (singletons) account for a large proportion of noisy reads. Because our emphasis was overall community trends (rather than exploration of the rare biosphere of corals), singleton sequences were removed. Representative sequences for each OTU were classified taxonomically according to the greengenes taxonomy version 13_866 using the Ribosomal Database Project (RDP) Classifier software v. 2.2. Sequence alignment and phylogenetic inference for the representative sequence of these OTU is described below in the context of β-diversity analysis.

We took additional steps to account for aspects of the data set unique to host-associated samples. Because coral mucus can contain some amounts of sloughed tissue, we tested whether coral mitochondria were present in any mucus samples. Similarly, because Symbiodinium and other photosynthetic microbial eukaryotes frequently inhabit coral mucus, chloroplast sequences are frequently observed in microbial diversity surveys of corals. As our interest was primarily in cellular bacteria and archaea rather than organelles, coral mitochondrial sequences and 16S sequences classified as chloroplasts with at least 70% confidence by the RDP classifier were removed before analysis.

Because we observed that many mitochondrial sequences were not efficiently identified by RDP, we also removed sequences with very high (1−e−50) sequence similarity and at least 90% sequence similarity to a reference coral mitochondrial sequence. A table of the coral mitochondrial sequences used is available as Supplementary Data 2 sheet b. These thresholds were selected with care to avoid indiscriminate removal of α-proteobacteria sharing evolutionary ancestry and 16S rRNA gene sequence similarity with mitochondria. The e-value for mitochondrial removal was selected by testing several e-values (1−e−10,1−e−30,1−e−50 and 1−e−100). For each e-value the best BLAST (Basic Local Alignment Search Tool) similarity against NCBI’s nr database was examined. We selected the highest (most lenient) e-value that removed mitochondria, but not related bacteria. At a 1−e−50 BLAST e-value threshold, the best BLAST hit of all removed 16S rRNA sequences was coral mitochondria, and this threshold was therefore selected for the screen. The final OTU table after quality control (QC) used for our analysis described in this manuscript can be found in Supplementary Data 7.

Analysis of microbial β-diversity

Microbial community β-diversity was calculated based on the weighted UniFrac distance matrix67. This is a phylogenetic measure of community similarity that takes into account organismal abundance and phylogeny67. Phylogenetic trees used to calculate this metric were constructed in QIIME 1.8 (ref. 65) through alignment of representative sequences of each OTU with PyNAST against the greengenes core set alignment66, and approximate maximum likelihood phylogenetic inference with FastTree. We considered the pool of distances between samples within each metadata category of interest (for example, algal competition or categories of temperature) using QIIME’s make_distance_boxplots.py. Significance was assessed by non-parametric t-tests, each with 1,000 Monte Carlo permutations (permutation is important in this instance to account for the non-independence of distances). The effect of this procedure is to ask whether different factors increase the dispersion of communities. PCoA plots of β-diversity were visualized in the Emperor software. When multiple categories (for example, different algal types) were tested for effects on β-diversity, the false discovery rate (FDR) for multiple comparisons was controlled at a threshold of q=0.05 using the Bejamini–Hochberg method.

Comparison of microbial and environmental distance matrices

Mantel tests test for correlation between two distance matrices. For example, a matrix of geographic distances for sample sites might be tested for correlation against a matrix of genetic distances. The partial Mantel test is an extension that tests for correlation between two distance matrices after accounting for the effects of a third, confounding, distance matrix. We used permutational Mantel tests to test whether between-sample variation in continuous environmental factors such as thermal stress, temperature or algal cover correlated with differences in the weighted UniFrac distance matrix67 between coral microbial communities. When data on hypothesized confounding factor was available, we used partial Mantel tests to test significance after accounting for the confounding parameter. For example, seasonal variation in algal cover might potentially confound the effects of temperature on microbial communities- partial Mantel tests were used to test for such effects. All Mantel and partial Mantel tests were performed in QIIME 1.8 (ref. 65) using the script compare_distance_matrices.py

Analysis of microbial community richness and evenness

Microbial community richness, evenness and β-diversity were calculated in QIIME v.1.8. Richness and evenness were calculated using the chao1 and equitability statistics, respectively, in QIIME’s alpha_diversity.py, collate_alpha.py and compare_alpha.ph scripts. In each case, the data were repeatedly (10 ×) rarified to 500 reads per sample. Values for chao1 and equitability were calculated for each rarified table, and averaged into a single value and compared across categories. Significance of results was assessed by FDR-corrected non-parametric t-tests with 1000 Monte Carlo permutations, using q<0.05 as the significance threshold.

Equitability was calculated as defined in the QIIME software package. Specifically, the Shannon entropy (using a base 2 logarithm) was divided by the base 2 logarithm of the number of observed OTUs. Thus, for a sample with n OTUs considered in turn (n 1 , n 2 , n 3 , …, n i ), the equitability/evenness was calculated as:

The effects of treatment, temperature, coral genus and individual coral head on richness and evenness were analysed using linear models in R (3.1.1). Seawater temperatures were binned into ‘low temperature’, ‘moderate temperature’ and ‘high temperature’ categories, with thresholds at<24 °C for ‘low temperature’ and >30 °C for ‘high temperature’. The high temperature threshold represents the MMM from the NOAA Pathfinder v5.2 SST data set45 (see above).

The specific form of these linear models in the R language was: lm(alpha_diversity ∼ treatment * seawater_ temperature_category * coral_genus+individual_coral_head). The significance of each effect was compared with ANOVA (Supplementary Data 3 sheet f).

Random forest analysis

Random forest analysis is a machine-learning method for supervised classification. We used random forest analysis (Supplementary Data 3 sheets i–j) to (i) determine the conditions under which Synechococcales would be displaced as the most abundant order in the coral microbiome, and (ii) generate a model predicting whether specific coral colonies would lose tissue using ecological (for example, long-term contact with algae from photo series) and microbiological data (the abundance of different bacterial orders).

For both analyses, random forest was conducted in QIIME v.1.8 with the script supervised_learning.py. In both models, the random forest was constructed with 1,000 decision trees and validated by 10-fold cross-validation. This analysis resulted in an overall prediction accuracy, a per-category prediction accuracy and a ranking of feature importance, defined as the amount of accuracy gained or lost when a specific feature was not included in the model. It is worth noting that we employed random forest analysis primarily to provide a summary of the overall predictive power of the data; the internal complexity of the algorithm does not allow inferences to be made about the relationships between variables, beyond their empirical influence on accuracy. For mechanistic insight we rely on other analyses presented here.

In the first model, we sought to predict the dominant (most abundant) order of bacteria in samples given only non-microbial data about the coral colony from which the sample was taken at that point in time. A total of 127 features were used in the prediction including: herbivore exclusion, nutrient addition, coral species (species were converted to binary variables), per cent cover of different macroalgae, turf algae and cyanobacterial mats from quadrat surveys, singly and grouped by functional type (for example, total macroalgae), plot and subplot number, measures of instantaneous and weekly average temperatures and salinities, season (quarter), and the simultaneous presence of parrotfish bites and nutrient enrichment. Given these features, we asked the model to predict the dominant order of microorganism in each DNA sample, with our main purpose being to determine whether the conditions under which Synechococcus was dominant were predictable. No microbiological data were included in the input feature table, so all predictions about the microbiology were based on externally observable features.

In the second model, we sought to test whether data about the competitive environment and microbial community composition of specific coral colonies could predict whether they would lose tissue over the 3-year study period. This model included features representing herbivore exclusion, nutrient pollution and coral genus; long-term competition with any algae or (as separate features) with turf algae, visible Cyanobacteria mats, Dictyota, Sargassum or Halimeda; the presence of bleaching, DSS or signs of parrotfish predation; the initial diameter of the coral colony in centimeters; and microbiological data (averaged across all samples from each individual) on the abundance of the phyla Proteobacteria, Cyanobacteria, Bacteroidetes, Firmicutes, Acidobacteria and Actionobacteria overall, and six specific microbial orders of interest due to their prominence in this data set and/or the literature13,41,44 (Synechococchales, Vibrionales, Rhodobacterales, Alteromonadales, Rickettsiales, Rhodospirallales and Pseudomonadales).

Predicting coral microbiome function

To predict potential functional consequences of observed changes in microbial taxonomy across treatments, functional profiles for each microbial sample were predicted using the PICRUSt tool40. This tool uses hidden state prediction, a form of evolutionary modelling closely related to ancestral state reconstruction, to put bounds on genomic copy numbers of each gene family in uncultivated environmental microorganisms, using their position in a reference bacterial phylogeny relative to all bacteria or archaea with sequenced genomes40,68. These per genome estimates are then combined, taking into account 16S rRNA copy number variation, to estimate a functional profile (or ‘virtual metagenome’) for a microbial community based on the information available in all sequenced genomes. The resulting metagenomes have been used in recent analyses of human and environmental microbiomes, and have in at least one study been shown to correlate with overall metabolite profiles (reviewed in ref. 68). The accuracy of the method depends on several factors40,68, but the availability of reasonably related reference genomes is generally the most important. Availability of such genomes can be summarized through a Nearest Sequenced Taxon Index (NSTI) score. This score is the average branch length between each member of the community (OTU) and its closest sequenced relative, weighted by abundance. In these data, NSTI scores ranged from 0.03 to 0.18, with a mean of 0.11. This was within the range of soil and mammal microbiomes that have previously been predicted with reasonable accuracy (as assessed using paired 16S rRNA/shotgun metagenomes from the same samples)40.

We used PICRUSt to predict functional profiles for all 16S rRNA samples, and summarized these profiles of predicted gene family abundance into Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways and functional categories. We then tested whether functional changes in coral microbiomes correlated with either temperature increase, extremes of temperature or the abundance of upright algal cover within a subplot (Supplementary Fig. 5a,b). To test the effects of extremes of temperature (either hot or cold) on coral microbiomes, the imputed abundance of KEGG functional categories (level 3) across samples was regressed against the squared deviation of temperature from 28 °C, which is the average temperature across samples, and also approximates the annual average temperature at the study site over 32 years (calculated from NOAA Pathfinder v5.245). The effects of upright algal cover or increasing temperature were tested by Spearman regression against KEGG functional category abundance. The FDR across KEGG categories was controlled in all cases at q<0.05.

Comparing PICRUSt results with previous experiments

To compare whether the functional categories predicted to change in this field study were broadly consistent with previous laboratory experiments, we compared PICRUSt data with a previous study that exposed corals to several stressors in aquaria and sequenced their metagenomes. KEGG pathways correlated with extremes of temperature (squared deviation from 28 °C, see ‘Predicting coral microbiome function’, above) in this study were compared with KEGG categories that increased or decreased in a previous experiment, in which Thurber et al.24, exposed the Pacific coral Porites compressa to thermal stress. Metagenomes from that study24 were annotated with KEGG pathways in MG-RAST using default parameters (BLAST e-value <10−5, 60% coverage and 15% alignable). We compared the 25 categories that increased by >1% with temperature stress to the present study.

Macroalgal contact as a driver of microbial β-diversity

We sought to test whether the increase in β-diversity with macroalgal contact that we saw in this data set is a common pattern in coral–algal competition. We reanalysed data from a previously published experiment19 that studied the effects of macroalgal contact on coral microbiomes of P. astreoides (Supplementary Fig. 7c). In that study, macroalgae were placed in direct contact with P. astreoides corals, and the microbial communities of the macroalgae, corals without algal competitors and corals with algal competitors assessed using terminal restriction fragment length polymorphism19. We translated the table of terminal restriction fragments from this study into a QIIME-compatible .biom format OTU table, and tested whether β-diversity (measured by Bray–Curtis divergence in terminal restriction fragment profiles) was greater in corals exposed to macroalgae than control corals or control algae using analyses described above (see Analysis of Microbial β-diversity).

Data availability

Processed microbial data are available as Supplementary Data 7 (OTU tables in BIOM format). Metadata for coral microbiome samples are available as Supplementary Data 2 sheets c,d, respectively. Short-read amplicon sequence data (as fasta and qual files) are deposited in the QIITA database, as study id: 10482. Other relevant data are available from the authors.