Multiple Trophic Tracer Analyses of Subarctic Rhodolith (Lithothamnion glaciale) Bed Trophodynamics Uncover Bottom-Up Forcing and Benthic-Pelagic Coupling

1 Introduction

Trophic ecology is the study of feeding relationships and energy transfers among organisms interacting in a community. In all ecosystems, energy is transferred through feeding from primary producers to primary and higher-order consumers. This transfer is often unclear in marine benthic ecosystems partly because of the broad diets of many species, a large detritus pool, and, sometimes, complex benthic-pelagic relationships, which can also vary seasonally (Kharlamenko et al., 2001; Pitt et al., 2009; Kelly and Scheibling, 2012). Trophic relationships can be studied with analysis of lipid classes, fatty acids (FA), and stable isotopes. Lipids are the densest form of energy (cal/g) in marine ecosystems, essential for structural integrity, storage, and signaling of molecules in cell membranes (Parrish et al., 2000; Parrish, 2009). Lipids and FAs can be used as biomarkers for food sources and nutrition because their composition can show input, cycling, and loss of material within food webs (Kelly and Scheibling, 2012). The intake, accumulation, and transferability of lipids and FAs make them an excellent tool to study trophic pathways (Richoux et al., 2005; Drazen et al., 2008b; Drazen et al., 2008a). Bulk stable isotope analysis is a useful approach to pair with FA analysis to help identify trophic relationships (Michener, 1994; Connelly et al., 2014), food and carbon sources (Carreón-Palau et al., 2013; Trueman et al., 2014), and food web structure (Grall et al., 2006; Linnebjerg et al., 2016).

Rhodoliths (free-living, non-geniculate red coralline algae growing as balls, branched twigs, or rosettes) often form dense aggregations, known as “rhodolith beds”, at depths of up to 150 m in tropical to polar seas (Foster, 2001; Foster et al., 2007). Rhodolith beds, along with seagrass meadows, kelp beds and forests, and mangrove forests, are one of the four major types of marine benthic primary producers (Foster, 2001; Foster et al., 2007). The relatively complex morphology of rhodoliths creates suitable habitats for attachment (Kamenos et al., 2004a; Steller and Cáceres-Martínez, 2009; Riosmena-Rodríguez and Medina-López, 2010), reproduction (Kamenos et al., 2004b; Steller and Cáceres-Martínez, 2009; Gagnon et al., 2012), and feeding (Steneck, 1986; Gagnon et al., 2012; Riosmena-Rodríguez et al., 2017) of highly diverse algal and faunal assemblages. The important contribution of rhodolith beds to marine biodiversity (Steller et al., 2003; Gagnon et al., 2012; Riosmena-Rodríguez et al., 2017) and global calcium carbonate (CaCO₃) production (Amado-Filho et al., 2012; Harvey et al., 2017; Teed et al., 2020) has, in part, triggered the recent increase in studies of factors and processes regulating their structure and function (Marrack, 1999; Hinojosa-Arango et al., 2009; Millar and Gagnon, 2018) and growth resilience to natural and anthropogenic stressors (Bélanger and Gagnon, 2020; Bélanger and Gagnon, 2021; Arnold et al., 2021)

Knowledge about trophodynamics in rhodolith beds is limited to only a couple of studies in northeastern Atlantic (Grall et al., 2006) and eastern Pacific (Gabara, 2014) systems, that together suggest suspended particulate organic matter (SPOM), sediment organic matter (SOM), and macroalgae are important components of rhodolith bed food webs. Both studies’ findings are based on use of bulk stable isotope analysis, in particular consideration of organisms’ carbon (δ¹³C) and nitrogen (δ¹⁵N) isotopic signatures (DeNiro and Epstein, 1978; DeNiro and Epstein, 1981; Minagawa and Wada, 1984) to identify primary producers (Peterson and Fry, 1987; Post, 2002; Bouillon et al., 2011) and trophic levels of consumers (Iken et al., 2001; Post, 2002). As noted by Newell et al. (1995) and Kelly and Scheibling (2012), benthic food webs with significant macroalgal and bacterial components are often too complex to be characterized uniquely with stable isotope analysis. This is likely the case for rhodolith beds, in particular those in seasonal seas, where phytoplankton blooms and growth of microalgal and bacterial films on the surface of benthic organisms occur seasonally. In such cases, FA and stable isotope analyses can, in principle, simultaneously help distinguish algal and bacterial inputs (Sargent et al., 1987; Kharlamenko et al., 2001; Kelly and Scheibling, 2012).

The present study aims to resolve a part of the trophic ecology of a large (>500 m²), well-studied, cold-water rhodolith (Lithothamnion glaciale) bed in southeastern Newfoundland (eastern Canada). Specifically, we combine stable isotope, lipid, and FA analyses to: (1) identify lipid compositions of organisms to better understand functional strategies in relation to environmental characteristics; (2) delineate trophic linkages among organisms to understand the nutritional value of their diets and the extent of benthic-pelagic coupling versus strictly benthic interactions; and (3) document specific challenges and requirements for future lipid- and stable isotope-based studies of feeding relationships in rhodolith beds. In doing so, we test the overall hypotheses that: (1) the lipid composition of rhodolith epifauna generally reflects the predominantly cold-ocean conditions of Newfoundland; and (2) the rhodolith bed food web is mainly controlled from the bottom up by planktivores and detritivores as reflected by high abundance of planktonic and bacterial biomarkers.

2 Materials and Methods

2.1 Study Site and Selection of Focal Species

Our study was carried out during the spring of 2017 in a rhodolith (Lithothamnion glaciale) bed, which extends ~5 to 30 m in depth along the coast of St. Philip’s, southeastern Newfoundland, Canada. Consistent with our broader objective of characterizing spatial and temporal variability in rhodolith bed trophodynamics, we chose to study the section of the bed fringing Broad Cove (47° 35′ 36.5″ N, 52° 53′ 31.0″ W; Figure 1A) because of presumed differences in marine environmental conditions in this area. Broad Cove is connected to a marina, which is the end point of a river originating from several large ponds in the centre of the northern Avalon Peninsula. Volumes of freshwater entering Broad Cove vary seasonally and are generally lowest in summer, when precipitations (ECCC, 2019), and hence river discharge, decrease (P. Gagnon, personal observations). Our study focuses on the southernmost of two sites, i.e. the “South” site (Figure 1A), which represents a relatively stable environment for the rhodoliths because of low hydrodynamic forces and sedimentation, and nearly unchanged rhodolith size structure over the past few years, as documented in the present and other studies of the bed (Figure 1B; Gagnon et al., 2012; Millar and Gagnon, 2018). We include the “North” site (Figure 1A) in a follow up study of trophic variability (Hacker Teper, 2022) because of the presumably greater environmental variability at this site resulting from greater proximity to freshwater input from the marina. The present study, therefore, focuses on trophic interactions at one site based on a single point in time, and serves as a baseline for expansion of spatial and temporal dimensions in a follow-up study (Hacker Teper, 2022).

Figure 1 (A) Location of the two study sites within the rhodolith (Lithothamnion glaciale) bed fringing St. Philip’s (southeastern Newfoundland) used to study rhodolith bed trophodynamics. The present study focuses on the food web at the “South” site in the spring of 2017. The “North” site is included in a follow up study (Hacker Teper, 2022) examining spatial and temporal variability in food web structure. Both sites are located at the periphery of Broad Cove, which receives seasonally variable volumes of freshwater from the adjacent marina and river to which it is connected (Image: Google Earth). (B) Section of the rhodolith bed at the South site at a depth of ~15 m. Rhodoliths are tightly aggregated, with very little to no epiphytes and a relatively high abundance of green sea urchins (S. droebachiensis) moving on the bed surface (the biggest urchins are ~6 cm in test diameter) (Image: Patrick Gagnon). (C) Staged photograph of rhodoliths (Lithothamnion glaciale) and associated macrofauna (visible: Ophiopholis aculeata, Strongylocentrotus droebachiensis, cryptic: Tonicella spp.) out of water from a laboratory bench at the Ocean Sciences Centre. (Image: Sean Hacker Teper).

Our food web analyses (described below) focused on the following six animal (1-6; Figures 1C, 2A−F), two macroalgal (7-8), and two environmental (9-10) components inside (1-6, 8, 10) or outside (7, 9) of the bed: (1) common sea star, Asterias rubens; (2) wrinkled rock-borer, Hiatella arctica [bivalve]; (3) juvenile Nereis spp. [polychaetes]; (4) daisy brittle star, Ophiopholis aculeata; (5) green sea urchin, Strongylocentrotus droebachiensis; (6) red molted chiton, Tonicella marmorea, and Atlantic red chiton, T. rubra, which were too difficult to distinguish morphologically, and hence were pooled to genus Tonicella spp.; (7) pieces of Laminaria digitata [kelp] from nearby kelp beds; (8) Lithothamnion glaciale [rhodoliths]; (9) seawater [containing seston] from a few meters above the rhodolith bed; and (10) sediment [containing infauna] underlying the rhodoliths. We chose these six animal species because they exhibit different dominant feeding strategies defined as per Macdonald et al. (2010), namely suspension/filter feeding (2, 4), grazing (5, 6), and predation (1, 3), but occur in sufficient abundance to quickly provide enough biological material for the analyses as per cryptofaunal diversity and abundance studies for this bed (Gagnon et al., 2012; Bélanger, 2020). Pieces of kelp (8), seawater from above the rhodolith bed (9), and sediment underlying the rhodolith bed (10) were sampled to explore possible benthic-pelagic coupling.

Figure 2 The six animal species included in the present study’s food web analyses; (A) common sea star, Asterias rubens; (B) wrinkled rock-borer, Hiatella arctica [bivalve]; (C) Nereis spp. [polychaete]; (D) daisy brittle star, Ophiopholis aculeata; (E) green sea urchin, Strongylocentrotus droebachiensis and (F) Atlantic red chiton, Tonicella rubra (Images: (A, B, D, F) - Sean Hacker Teper; (C) https://www.enasco.com/p/Sandworm-Clam-Worm-Nereis%2C-Preserved%2BLS01292; (E) https://pugetsoundsealife.sseacenter.org/pugetsoundsealife.com/puget_sound_sea_life/Green_Sea_Urchin.html.

2.2 Timing of Sampling

To increase the likelihood of detecting benthic-pelagic coupling, we sampled the rhodolith bed during the annual spring phytoplankton bloom in southeastern Newfoundland, when diatom abundance in the water column was predictably highest (Budge and Parrish, 1998; Parrish et al., 2005). We collected fluorescence data at the study site with a CTD (conductivity, temperature, and depth) profiler equipped with PAR (photosynthetically active radiation) and fluorescence sensors to monitor the progression and confirm the occurrence of the bloom (Appendix A). The bloom began in the last few days of March 2017 and continued until at least 23 April, 2017 (Figure A.1A), when we sampled the rhodolith community and collected rhodoliths for food web analyses.

2.3 Rhodolith Community

To broadly characterize the rhodolith community at our collection/study site, scuba divers hand collected, on 23 April, 2017, all the rhodoliths from one 30 x 30 cm quadrat placed every 5 m along a 30-m long transect at a depth of ~15 m at the South site (for a total of seven quadrats sampled). We chose to sample this section of the bed because of the relatively homogenous distribution and high abundance of rhodoliths there (Gagnon et al., 2012; Millar and Gagnon, 2018); Figure 1B). Rhodoliths and their associated cryptofauna were deposited in labelled, sealable plastic bags, with several bags per quadrat. Bags were sealed under water, placed in mesh collection bags, and lifted to a boat where they were stored in plastic bins filled with seawater collected at the site. They were transported to the Ocean Sciences Centre (OSC) of Memorial University of Newfoundland (MUN) within 4 h of collection and placed in large (320-L) holding tanks supplied with running seawater pumped in from the adjacent Logy Bay.

Lengths of the longest, intermediate, and shortest axes, as well as gross weight of the 247 rhodoliths collected were measured with calipers (precision of 0.1 mm) and a balance (precision of 0.1 g; PB3002-S/FACT; Mettler Toldeo). Each rhodolith’s linear dimensions and number of rhodoliths in each of the seven quadrats were subsequently used to calculate each rhodolith’s sphericity (as per Graham and Midgley [2000] and Sneed and Folk [1958]) or to estimate rhodolith abundance (density) in the bed. We plotted rhodolith shape with the spreadsheet TRIPLOT (https://www.lboro.ac.uk/microsites/research/phys-geog/tri-plot/tri-plot_v1-4-2.xls) (Sneed and Folk, 1958; Graham and Midgley, 2000) as described by Gagnon et al. (2012). All cryptofauna on the external surface and inside of each rhodolith were extracted with tweezers and forceps, breaking rhodoliths in pieces with a screwdriver and a hammer to extract hidden specimens when needed. We placed organisms in labeled specimen cups filled with a 4% formalin solution prior to permanent changeover after three days into a 70% ethanol solution. Over the following few weeks, preserved cryptofauna were identified and counted with a stereomicroscope (DMW-143-N2GG; Motic) at 10 or 20X magnification and weighed with a balance (precision of 0.1 g; PB3002-S/FACT; Mettler Toldeo). Epiphytes and encrusting invertebrates such as bryozoans and sponges were present in trace amounts, and hence excluded from the analysis. For each quadrat, we subsequently subtracted total cryptofaunal weight from gross rhodolith weight to obtain net rhodolith weight, which we used to calculate rhodolith abundance (biomass) in the bed.

2.4 Collection and Preparation of Samples for Food Web Analyses

To limit influences of manipulation of rhodoliths and their cryptofaunal content on data quality, we used a different group of rhodoliths for food web analyses than those sampled to characterize the rhodolith bed community. On 23 April, 2017, divers hand collected ~150 live rhodoliths measuring 8 to 10 cm along the longest axis from the same area where we sampled the associated bed community. Broken rhodoliths and rhodoliths partially buried in sediment or with blackened or whitened tissue on their surface (indicative of stress or necrosis) were not collected because of potential influences on the abundance, diversity, and chemical composition of resident cryptofauna. Preliminary analysis indicated that common sea star and green sea urchin biomasses within the rhodoliths were too low to provide the minimum amount of tissues required to run the lipid and stable isotope analyses (see below). We resolved this requirement by collecting an additional ~10 small individuals (1-2 cm in diameter) from atop the bed for each species. The fronds of three, ~1 m long kelp (L. digitata) growing on rocks at a depth of ~2 m near the rhodolith bed were also hand collected. Rhodoliths and kelp fronds were placed in plastic bags sealed underwater. We collected seawater a few centimeters above the rhodolith bed with two, 12-L Niskin bottles that we deployed gently to prevent resuspension of sediment from the bed. Water from the bottles was transferred (on the boat) to plastic containers pre-rinsed with distilled water, from which a total of 17 L of seawater was subsequently taken to meet the requirements of the various analyses (see below). Three sediment samples were also scooped from the top (~10 cm) layer of muddy sediment underneath rhodoliths with 15-mL centrifuge tubes.

All rhodoliths and their cryptofaunal content, as well as kelp fronds, were transported to the OSC as described above (see section 2.3). At the OSC we transferred all the rhodoliths in their sealed plastic bags to large (320-L) holding tanks supplied with running seawater to keep water temperature in the bags naturally low (~0.5°C), while retaining all cryptofauna. Bags were kept sealed in the tanks for 24 h to (1) facilitate stomach emptying of focal species; (2) prevent hypoxia and degradation of biological tissues; and (3) avoid mixing water from St. Philip’s [in the bags] with water from Logy Bay [in the holding tanks]. This procedure reduced the likelihood of contaminating the original lipid class, FA, and isotopic signatures of the focal species. After this 24-h resting phase, we removed bags from the tanks and processed them one by one with pre-cleaned/sterilized tools and glassware manipulated with frequently changed nitrile gloves, again to avoid contamination of samples. Test tubes and scintillation vials were heated in an oven for 8 h at 425°C to remove lipid material, then labeled and weighed. Test tube caps and the tools used to break rhodoliths (hammer, screwdriver, mortar, and pestle) and extract and separate cryptofauna (tweezers, forceps, filters) were lipid-cleaned with three rinses each of methanol and chloroform to remove any residual lipids. We rinsed collection and storage items (sealable plastic bags and containers, centrifuge tubes) with distilled water.

We extracted the six focal animal species from the surface of rhodoliths with tweezers and forceps, breaking rhodoliths in pieces with a screwdriver and a hammer to extract cryptic specimens when needed. For lipid class and FA analyses, we obtained three replicates of 0.5 to 1.5 g of tissue each (wet weight) for each species from one or several individuals (i.e. pooling tissue as needed). We included whole individuals (i.e. shells, exoskeletons, and internal organs) in the wet weight of lipid class and FA analysis samples. Each replicate was placed in a 40-mL test tube (one replicate per tube) and stored on ice in a cooler until we had collected all replicates. Test tubes were then filled with 8 mL of chloroform, flushed under a gentle stream of nitrogen, capped, sealed with Teflon tape, stored in freezers at -20°C, and their content analyzed within two months. For stable isotope analysis, we collected ~5 g of tissue, excluding shells and exoskeletons, into one 20-mL scintillation vial for each species to be split into three replicates once dry (see section 2.7). Due to the randomness of availability of material of the small, pooled organisms and because lipid extraction may affect δ¹⁵N data (Post et al., 2007), we did not use the same organism material between lipid-extraction and stable isotope analysis (see section 2.9 for statistical implications). After oven drying for 24 h at 60°C, at least 1 to 1.5 mg of tissue remained for each replicate. Each scintillation vial was capped with tin foil and a cap, stored in freezers at -80°C, and their content analyzed within 12 months.

Upon arrival at the OSC, kelp fronds and sediment samples were immediately stored in their individual, sealed plastic bags (kelp) or centrifuge tubes (sediment) in freezers at -80°C. We used the same procedures as above to prepare rhodolith, kelp, and sediment samples for lipid class, FA, and stable isotope analyses, with the following modifications. Rhodoliths were first gently scrubbed by hand to remove epibionts, broken off with a screwdriver and a hammer, and ground into a powder with a mortar and a pestle. Rhodolith powder was analyzed; within two (lipid class and FA analyses) or 12 (stable isotopes) months. Kelp fronds (blades and stipes) were gently scrubbed by hand to remove epibionts. Each replicate of sediment for the lipid class and FA analyses weighed 6 to 8 g (wet weight). We also analyzed kelp and sediment samples within two (lipid class and FA analyses) or 12 (stable isotopes) months.

We processed seawater samples immediately upon arrival at the OSC. The two plastic containers holding the seawater were shaken to re-suspend any settled materials. We transferred this water and its content with 250- to 500-mL graduated cylinders to a mechanized filtration system, which suctioned water with an aspirator through a 47-mm diameter GF/C filter (Whatman; General Electric) at the bottom of a Büchner funnel. During suction, we washed the contents of graduated cylinders onto the filters with filtered seawater to transfer all lipid material. Visual inspection of the first filter indicated an acceptable accumulation of suspended materials upon completing the filtration of the first 3 L of seawater. Given the total volume of water available (~17 L), we created three replicates for the lipid class and FA analyses; one per each 3 L of filtered seawater. We also filtered the remaining water, 8 L, yielding three replicates for the stable isotope analysis; one per each 2.5 L of filtered seawater. Each filter used was rolled with tweezers and placed in a 40-mL test tube (one filter per tube) for the lipid class and FA analyses, or in a 20-mL scintillation vial (one filter per vial) for the stable isotope analysis. Test tubes were filled with 8 mL of chloroform, flushed under a gentle stream of nitrogen, capped, sealed with Teflon tape, stored in freezers at -20°C, and their content analyzed within two months. Vials were capped with tin foil and a cap, stored in freezers at -80°C, and their content analyzed within 12 months.

2.5 Extraction and Characterization of Lipid Classes

Extraction of lipids followed protocols by Folch et al. (1957) with modifications by Parrish (1999). We took samples in their test tubes out of the freezers and handled 8 to 16 at a time. Each tube was held in ice and contents were immediately ground to a pulp with a metal-ended rod washed into the tube with ~1 mL of chloroform:methanol (2:1) and 0.5 mL of chloroform-extracted water. We sonicated the tube for 4 min and centrifuged it for 3 min at 3000 rpm. The resulting organic layer at the bottom of the tube was completely removed and transferred to a 15-mL vial with a double pipetting technique to bypass the upper aqueous layer and transfer only the lower aqueous layer. After three repetitions, we washed both pipettes into the 15-mL vial with 3 mL of chloroform. All vials were flushed under a gentle stream of nitrogen, capped, sealed with Teflon tape, and stored in freezers at -20°C. We subsequently transferred each vial’s content to a 2-mL, lipid-clean vial with three or more rinses of 0.5 mL of chloroform, until the chloroform in the 15-mL tube remained transparent. Each 2-mL vial was flushed under a gentle stream of nitrogen, capped, sealed with Teflon tape, and stored in freezers at -20°C.

We used thin-layer chromatography with flame ionization detection (TLC-FID) to characterize lipid classes (Parrish, 1987). Separation of lipid classes followed a 3-step development method in which four solvent solutions of different polarities were used to obtain three chromatograms per rod (Parrish, 1987). In a first step, we created calibration curves (against which samples were compared) with a 9-component standard (nonadecane, hydrocarbon; cholesteryl palmitate, steryl ester; 3-hexdecanone, ketone; tripalmitin, triacylglycerol; palmitic acid, free fatty acid; cetyl alcohol, alcohol; cholesterol, sterol; monopalmitoyl glycerol, acetone mobile polar lipid; phosphatidylcoline dipalmitoyl, phospholipid; Sigma Chemicals). Prior to applying samples, silicic acid-coated quartz rods (Chromarods, Type SV; Iatron Laboratories Inc.) were blank scanned (the process of cleaning and activating rods by burning off any residual lipids from previous samples) three times in an Iatroscan TLC-FID system (Mark VI; Iatron Laboratories Inc.). We operated the Iatroscan with a hydrogen flow between 195 and 199 mL min^-1 and an air flow of 2 L min^-1 at a room temperature of ~20°C. After sample application (spotting) and before each development, the rods were dried and conditioned for 5 min in a constant humidity chamber (30%). We spotted samples individually onto one of 10 rods in each of 2 racks (for a total of 20 rods) with a 25-μL Hamilton syringe (Hamilton Co.). Depending on the lipid concentration of each sample (estimated by colour and confirmed by trial and error), we spotted 0.5 to 10 μL was spotted at the origin of each rod. Both racks were then dipped twice into a 100% acetone solution and were removed each time when acetone reached the spotted samples. We then double-developed rods in a hexane:diethyl ether:formic acid (98.95:1:0.05) solution, for 25 min, and then for another 20 min to separate hydrocarbons (HC), steryl esters (SE), and ethyl and methyl ketones (KET). We scanned each rod in the Iatroscan beyond the KET peak to obtain the first chromatogram showing lipid quantities of HC, SE, and KET in each sample.

In a second step, we developed rods for 40 min in a hexane:diethyl ether:formic acid (79:20:1) solution to separate diacyl glyceryl ethers, triacylglycerols (TAG), free fatty acids (FFA), alcohols (ALC), sterols (ST), and diacylglycerols (DG). We scanned each rod beyond the DG peak to obtain the second chromatogram showing lipid quantities of TAG, FFA, ALC, ST, and DG in each sample. In a third step, we double-developed rods twice, first in a 100% acetone solution for 15 min each, then in a chloroform:methanol:chloroform-extracted water (5:4:1) solution for 10 min each to separate the most polar lipid classes, acetone-mobile polar lipids (AMPL), and phospholipids (PL). Finally, each rod was completely scanned to obtain the third and final chromatogram showing lipid quantities of AMPL and PL.

2.6 Preparation and Characterization of Fatty Acid Methyl Esters (FAME)

We prepared fatty acid methyl esters (FAME) of lipids by transesterification of extracts following a modified procedure described by Christie (1982) and Hamilton (1992). We took extracted lipid samples in their 2-mL vials from the freezer and held them in ice. Depending on each sample’s lipid concentration (estimated by colour and confirmed by trial and error), we transferred 20 to 1000 µL of lipid extract with a 20 to 100 µL Drummond microdispenser (Drummond Scientific) into a lipid-clean, 15-mL vial. We evaporated the transferred extract to dryness under a stream of nitrogen gas, then we added 1.5 mL of dichloromethane and 3 mL of prepared Hilditch reagent (an alkylation derivatization reagent, 1.5 H₂SO₄: 98.5 MeOH). Each vial was vortexed, sonicated for 4 min, flushed with nitrogen, heated at 100°C for one hour, cooled to room temperature, filled with 0.5 mL of a supersaturated sodium bicarbonate solution and 1.5 mL of hexane, and agitated vigorously. This process created an organic layer containing fatty acids transesterified to fatty acid methyl esters (FAME), which we transferred with a pipette to a lipid-clean, 2-mL vial, and evaporated to dryness under a stream of nitrogen. Each 2-mL vial was then filled with ~0.5 mL of hexane, flushed with nitrogen, capped, sealed with Teflon tape, and stored in freezers at -20°C.

We used gas chromatography and flame ionization detection (GC-FID) to measure fatty acids (Ackman, 1986; Christie, 1989; Budge et al., 2006). FAME samples were analyzed in a HP 6890 GC equipped with an Agilent 7683 autosampler. The GC column was a 30-m long ZB wax+ (Phenomenex) with an internal diameter of 0.32 mm. We introduced each FAME sample individually into a heated injector at 150°C. A stream of hydrogen then carried the sample at a rate of 2 mL min^-1 through the GC column, which retained FAME depending on structure. This selective retention resulted in the detection of individually eluted FAME by the FID. The GC column temperature started at 65°C for 30 s before it was increased to 195°C at a rate of 40°C min^-1, held there for 15 min, then increased again to 220°C at a rate of 2°C min^-1, and held there for 45 s. The initial injector temperature of 150°C increased to a final temperature of 250°C at a rate of 120°C min^-1. The detector (FID) temperature remained constant at 260°C. We identified FA peaks by comparing retention times in the GC column with those from various standards: FAME mix (47885-U; Supelco), Bacterial acid methyl ester mix (47080-U; Supelco); PUFA 1 (47033; Supelco); and PUFA 3 (47085-U; Supelco). We integrated chromatogram peaks with Varian Galaxie Chromatography Data System V1.9.3.2 software.

2.7 Stable Isotope Preparation and Analysis

We removed samples in scintillation vials from the freezer and thawed and oven dried them at 60°C for 24 h. Each sample was split equally by weight into three new, 20-mL scintillation vials. Each triplicate was ground with a mortar and pestle and rinsed three times with acetone. We then split each ground sample into two equal parts: one for analysis of carbon content, the other for nitrogen. Vials with samples for nitrogen analysis were each rinsed three times with distilled water to remove salt because halides interfere with the Elemental Analyzer (described below). We acidified samples for carbon analysis to remove inorganic carbon by adding drops of 1 M HCl to the vials until no bubbles formed. Vials were left with lids off in a fume hood overnight, rinsed three times with distilled water to remove both HCl and halides, dried again at 60°C for 24 h, transferred to desiccators, and subsequently taken to The Earth Resources Research and Analysis (TERRA) facility at MUN for analysis. We used tools and tin capsules pre-cleaned with acetone to extract, weigh, and hold the following quantities of solidified samples from each vial: 1 to 1.5 mg for animals; 4 to 4.5 mg for rhodolith and kelp; 14 to 15 mg for sediment; 5 to 7 mg of filter for seawater δ¹³C and carbon content; and 8 to 10 mg of filter for seawater δ¹⁵N and nitrogen content. We held tin capsules and their content in a desiccator and processed them within a month as per the following procedure.

We analyzed bulk stable isotopes and complementary C and N elemental proportions (%) in an Elemental Analyzer (EA) system (NA1500; Carlo-Erba) consisting of an autosampler, an oxidation reactor (oven), a reduction reactor, a water trap, a gas chromatography (GC) column, and a thermal conductivity meter (TCD). The entire EA was flushed continuously with helium at a rate of 90 to 110 mL min^-1. We dropped each tin capsule and its content individually onto the oxidation reactor at a temperature of 1050°C, with simultaneous injection of oxygen and quick flushing with He. This sequence triggered a flash combustion at 1800°C between the tin capsule and oxygen, creating combustion gases that were pushed through an oxidation catalyst (chromium trioxide, CrO₃) to ensure complete oxidation of the sample and silvered cobaltous/cobaltic oxide, which removes halides and SO₂. The resulting gas mixture passed through the reduction reactor (reduced copper) at 650°C, which reduces nitrogen oxides to nitrogen gas and absorbs oxygen. The gases then passed through a magnesium perchlorate Mg(ClO₄)₂ water trap, after which the remaining gases (N₂, CO₂) entered a 3-m GC column (QS 50/80; Poropak) at 40 to 100°C. The individual gases were on the GC column. Upon reaching the TCD, they were detected as separate gas peaks; first N₂, then CO₂. From the TCD, He carried the gases to a ConFloIII interface (Finnigan, Thermo Electron Corporation), which has split tubes, open to the atmosphere, which allow a portion of the He and combustion gases to enter directly into the ion source of the mass spectrometer (MS) (DeltaVPlus; Thermo Scientific). During operation, He from the EA flowed continuously into the MS. Internal and external reference material was used to calibrate MS data. We used EDTA #2 and D-Fructose for carbon isotope calibration, and IAEA-N-1 ((NH₄)₂SO₄) and IAEA-N-2 ((NH₄)₂SO₄) for nitrogen isotope calibration. We used NBS-18 (CaCO₃), B2150 (high organic sediment), B2151 (high organic sediment), and B2105 (Cystine) to aid data interpretation of carbon isotope analyses, and sorghum flour, B2153 (low organic soil), USGS-25 ((NH₄)₂SO₄), USGS-26 ((NH₄)₂SO₄), sulfanilamide, and BBOT to aid data interpretation of nitrogen isotope analyses. L-glutamic acid and B2155 (protein) were used for both carbon and nitrogen elemental calibration.

2.8 Trophic Magnification of Fatty Acids

Stable isotope ratios are expressed in the conventional (δ) notation as parts per thousand (‰) as per the equation of Minagawa and Wada (1984):

where R_sample and R_standard are the ratios of ¹³C/¹²C or ¹⁵N/¹⁴N of a given sample and corresponding standard, respectively. Results are reported relative to atmospheric N₂ for nitrogen stable isotopes, and Vienna PeeDee Belemnite (VPDB) for carbon stable isotopes. Species trophic position (TP_consumer) was calculated with the equation used by Gale et al. (2013), and developed originally by Cabana and Rasmussen (1996):

where δ¹⁵N_consumer is the mean stable N isotope ratio of each species, and Δ¹⁵N is the fractionation factor which, to be consistent with rhodolith food web studies, is 3.4‰ (Grall et al., 2006). δ¹⁵N_base and TP_base represent the nitrogen stable isotope composition and trophic positions from the base of the food web, respectively. We then calculated a trophic magnification factor (TMF) for fatty acids (FA) correlated with δ¹⁵N. This factor quantitatively represents the biomagnification of compounds along a food web (Borgå et al., 2012; Connelly et al., 2014). Compound concentrations often change across trophic levels, thus the equation:

or

and, therefore,

where m and b are the slope and intercept of the linear relationship (which strength was determined with conventional Spearman Rank-Order correlation tests) between log_e FA% and trophic position (TP), respectively. Positive values of m and TMF imply biomagnification throughout the food web, whereas negative values denote proportional depletion. To this extent, we did not use TMF as a tool to see changing concentrations of FAs through direct predator-prey relationships among organisms. Instead, we used TMF as a method to assess biomagnification and depletion of FAs by comparing FAs to one another and to help confirm the presence of key identified biomarkers throughout the food web (Connelly et al., 2014).

2.9 Statistical Analysis

2.9.1 Lipid Classes

We used a one-way permutational MANOVA (PERMANOVA) (Euclidean distance matrices with 9999 permutations) with the factor Component (nine of the 10 components of the food web studied [six animal species, kelp, sediment, and seawater]) to examine differences in proportions of lipid classes among samples (N=25, accounting for accidental loss of two samples during the analyses). We excluded the rhodolith component because of insufficient rhodolith tissue for lipid extraction. To limit extraneous data variability while focusing on the most significant lipid classes, our analysis included only lipid classes present in over 50% of the samples. Sample sizes for each food web component (N=2 or 3) was too low to examine differences among specific components. Consequently, for comparison purposes only, we pooled the data into the five following functional groups reflecting the three dominant feeding strategies of the six animal species, one macroalgal species, and two environmental components: (1) suspension/filter feeders [two species]; (2) grazers [two species]; (3) predators [two species]; (4) kelp; and (5) seawater/sediment [samples combined because of expected benthic-pelagic coupling and to achieve a sufficient sample size for statistical analysis] (see section 2.1 for details). We then ran a one-way PERMANOVA with the factor Functional Group (the five groups explained above). This approach yielded statistically reliable comparisons, except with kelp, for which sample size was too low and could not be pooled with any of the four other groups because of its unique nature. We therefore do not present comparisons with kelp. We examined relationships between total lipid and each of the major lipid classes with conventional Spearman Rank-Order Correlation tests (Zar, 1999).

2.9.2 Fatty Acids

To examine differences in the proportions of FAs among samples (N=25), we used the same statistical approach (two one-way PERMANOVAs; one with Component as factor followed by one with Functional Group as factor) with the same data exclusion and grouping as for the lipid classes analysis. We then used a one-way SIMPER analysis (run on untransformed data with a Bray-Curtis similarity matrix) with the factor Component (nine of the 10 components of the food web studied [six animal species, kelp, sediment, and seawater], to identify potential food sources and the main FAs contributing to the lipid composition of each component (Kelly and Scheibling, 2012; Gabara, 2014). To limit extraneous data variability while focusing on the most significant FAs, we included only FAs contributing to over 70% of the similarities in the SIMPER analysis. We used a follow-up principal coordinates analysis (PCO; also run on untransformed data with a Bray-Curtis similarity matrix) with the factor Component (same as above), mainly for visualization of the feeding relationships among specific groups of organisms (Guest et al., 2008; Drazen et al., 2009). To increase clarity on the PCO, we plotted only samples with a Pearson coefficient of correlation >65% (plus DHA; 22:6ω3, because of its importance as an essential FA).

2.9.3 Stable Isotopes

We examined differences in carbon (δ¹³C) and nitrogen (δ¹⁵N) isotope ratios with two one-way PERMANOVAs (one for each type of ratio; both types based on Euclidean distance matrices with 9999 permutations) with Component as factor. The δ¹³C isotope ratio analysis (N=27, accounting for accidental loss of three samples during the analyses) included all taxa (six animal species, kelp, and rhodoliths) and both environmental components (seawater and sediment). We included both environmental components and all except one taxon (Nereis spp., for which we lacked sufficient amounts of tissues for quantification of stable isotope ratios) in the δ¹⁵N isotope ratio analysis (N=26, accounting for accidental loss of one sample during the analyses). Due to our using separate samples for lipid-extract and stable isotope analyses, we used the averages of stable isotope results of individual components to make comparisons to lipid and FA data. To group and map, in the form of a dendrogram, statistically different components of the food web, we also carried out a cluster analysis using “Group Average” clustering on δ¹³C and δ¹⁵N isotope ratios simultaneously, and complementary SIMPROF test (Euclidian distance matrix with 9999 permutations) (N=23 because of a few unmatched pairs of δ¹³C and δ¹⁵N ratios) (Grall et al., 2006; Gabara, 2014). Four main isotopic groups emerged from the SIMPROF test. We therefore ran follow-up one-way PERMANOVAs (one with both isotopic ratios combined, followed by one for each type of isotopic ratio) and complementary one-way ANOVAs and post-hoc tests with the factor Group to identify differences among these four main trophic groups.

2.9.4 General Aspects of Statistical Tests

In all PERMANOVAs, data were untransformed and computed on Bray-Curtis similarity or Euclidian distance matrices (9999 permutations) to meet the assumptions of multivariate normal (Gaussian) distribution and homogeneity of the covariance matrices (Budge et al., 2006; Clarke et al., 2006; Hair et al., 2006). All FA multivariate data were computed using Bray-Curtis similarity matrices, while lipid and stable isotope multivariate data were computed using Euclidean distance matrices due to its better ability to handle missing data. Contrary to the recommendation from Kelly and Scheibling (2012), we used untransformed data because the dispersion of variance was equivalent to transformed data and it avoided artificial weighting of FAs with smaller proportions on our results compared to transformed data (Carreón-Palau et al., 2017). We used PERMDISP (9999 permutations) to inform our decision (p=0.391); we tested for homogeneity of multivariate variances and confirmed all variances were homogenous. We used PCO (principal coordinates analysis) instead of PCA (principal components analysis) to more efficiently account for missing data (Rohlf, 1972). In all ANOVAs, we verified homogeneity of variance and normality of residuals by examining the distribution of the residuals and the normal probability plot of the residuals, respectively (Snedecor and Cochran, 1994). We used a significance level of 0.05 in all analyses and report all means with standard deviation (mean ± SD) unless stated otherwise. We used standard error where applicable for consistency with corresponding literature (Gagnon et al., 2012; Connelly et al., 2014; Parzanini et al., 2018). We used PRIMER v7 with PERMANOVA+ for multivariate statistical analyses, Minitab 18 for univariate statistical analyses, and Microsoft Excel for descriptive statistics.

3 Results

3.1 Rhodolith Community

Rhodolith biomass at the study site averaged 19.5 ± 0.9 (SE) kg m^-2. The 247 rhodoliths sampled varied in size from 11.3 to 65.6 mm, and 24.0 to 116.2 mm along the shortest and longest axes, respectively. Those rhodoliths were predominantly spheroidal and compact (~61%), but otherwise platy (~12%), bladed (~19%), or elongate (~8%) Figures 3A, B). Total cryptofaunal biomass averaged 34.5 ± 4.3 g kg^-1 rhodoliths. The 1191 animals extracted from the rhodoliths belonged to at least 21 species under six phyla, with echinoderms (452.7 ± 47.0 individuals kg^-1 rhodoliths) and molluscs (427.6 ± 39.6 individuals kg^-1 rhodoliths) as the two numerically dominant groups (Table 1). Species included in the biochemical analyses were particularly abundant, including Ophiopholis aculeata (336.7 ± 30.8 individuals kg^-1 rhodoliths), Tonicella spp. (191.6 ± 27.2 individuals kg^-1 rhodoliths) and Hiatella arctica (152.8 ± 23.2 individuals kg^-1 rhodoliths). A few species not included in the analyses were also relatively abundant, including the brittle star, Ophiura robusta (72.7 ± 13.8 individuals kg^-1 rhodoliths), the caridean shrimp, Pandalus borealis (36.9 ± 7.5 individuals kg^-1 rhodoliths), and the polychaete, Potamilla reniformis (30.0 ± 10.3 individuals kg^-1 rhodoliths).

Figure 3 (A) Ternary diagram of rhodolith [Lithothamnion glaciale] shape relative to purely spheroidal, discoidal, and ellipsoidal rhodoliths [N=247; one solid circle per rhodolith]. Rhodoliths were collected in April 2017 at the South site (see Figure 1A). The position of each rhodolith in the diagram is determined by its sphericity, calculated from the length of its longest [L], intermediate [I], and shortest [S] axes. (B) Corresponding proportion of rhodoliths under each of 10 finer shape categories as defined by Sneed and Folk (1958).

Table 1 Taxonomical breakdown and abundance of invertebrate cryptofauna associated with rhodoliths (Lithothamnion glaciale) collected in April 2017 at the South site (see Figure 1A).

3.2 Total Lipid Content and Lipid Classes

Of the nine food web components included in the lipid analysis, the three echinoderm species exhibited the highest concentrations of total lipids, ranging from 8.5 ± 2.2 mg g^-1 (ww) in A. rubens, to 13.2 ± 3.7 mg g^-1 in O. aculeata (Table 2). The molluscs H. arctica and Tonicella spp. had, respectively, similarly high and slightly lower concentrations of total lipids, with 6.7 ± 1.6 mg g^-1 in Tonicella spp. (Table 2). The polychaete Nereis spp. exhibited the lowest concentration among animals, with 6.0 ± 1.9 mg g^-1. Kelp (L. digitata), seston (from seawater samples), and sediment had significantly lower total lipid concentrations than O. aculeata (Tukey HSD, p<0.01 in all cases) (Table 2). The nine food web components contained nine lipid classes (PL, TAG, FFA, ST, AMPL, HC, SE, KET, and ALC), with six (PL, TAG, FFA, ST, AMPL, and HC) present in >50% of all samples (Table 2). PL was the dominant lipid class in every component, with a proportional contribution to total lipid concentration of 48% in O. aculeata to 76% in A. rubens (Table 2). Animal species contained the highest proportion of TAG, ST, and AMPL, with 34% in H. arctica, 19% in Nereis spp., and 16% in S. droebachiensis, respectively. FFA was highest (31%) in seston, lowest (2%) in A. rubens, and not detected (0%) in S. droebachiensis (Table 2). Kelp and seston contained no measurable TAG. Seston lipids were largely PL (50%), FFA (31%), ST (11%), and AMPL (12%). Sediment was dominated by PL (53%), followed by TAG (16%) and AMPL (13%). All components contained PL, ST, and AMPL. Overall, lipid class composition differed significantly among the nine food web components (PERMANOVA; Pseudo-F_{8,24 =} 5.732, P (perm)=0.0003) and five functional groups (suspension/filter feeder, grazer, predator, kelp, and seawater/sediment) studied (PERMANOVA; Pseudo-F_4,24= 3.504, P (perm)=0.0059). Specifically, lipid class composition differed between grazers and predators (t=2.485, p=0.018), and between suspension/filter feeders and predators (t=4.450, p=0.003). Total lipid concentration correlated with ST proportion only (r=-0.469, p=0.018, N=25).

Table 2 Sample size (N), mean wet weight, mean total lipid, and mean proportion (%) of the six dominant lipid classes (PL, phospholipid; TAG, triacylglycerol; FFA, free fatty acid; ST, sterol; AMPL, acetone mobile polar lipid; and HC, hydrocarbon) in the six animal species (common sea star, Asterias rubens; wrinkled rock-borer, Hiatella arctica [bivalve]; Nereis spp. [polychaetes]; daisy brittle star, Ophiopholis aculeata; green sea urchin, Strongylocentrotus droebachiensis; and chitons, Tonicella spp.), two macroalgal species (Laminaria digitata [kelp] and Lithothamnion glaciale [rhodoliths]), and two environmental components (seawater and sediment) sampled inside (I) or outside (O) of the South site (see Figure 1A).

3.3 Fatty Acid Profiles

The nine food web components included in the FA analysis contained 63 FA, with 43 present in >50% of all samples. Each component exhibited a distinct set of dominant FAs and biomarkers (Table S.1 and Figure 4). With a proportional contribution to FA profiles ranging from 19% in H. arctica to 32% in A. rubens, EPA (eicosapentaenoic acid, 20:5ω3; a typical diatom biomarker) was the dominant FA within each of the six animal species sampled (Table S.1). ARA (arachidonic acid, 20:4ω6; a kelp and amphipod biomarker), palmitoleic acid (16:1ω7; a diatom biomarker), and palmitic acid (16:0; a flagellate, bacteria, and marine vegetation biomarker) were the next most prominent FA among the animals, with a contribution between 11% and 20% (Table S.1). Kelp contained mainly palmitic acid (18%), eicosatetraenoic acid (20:4ω3; 14%) and oleic acid (18:1ω9, 13%; a crustacean, detritus, dinoflagellate, and brown seaweed biomarker). Seston FA were largely oleic acid (33%), palmitic acid (22%), and stearic acid (18:0, 18%; a detritus biomarker). Sediment was dominated by palmitoleic acid (21%), palmitic acid (15%), and vaccenic acid (18:1ω7, 9%; an aerobic microorganism, bacteria, and vegetation biomarker) (Table S.1). Overall, FA composition differed significantly among the nine food web components (PERMANOVA, Pseudo-F_8,24 = 26.278, P (perm)<0.001) and five functional groups studied (PERMANOVA, Pseudo-F_4,24 = 7.6664, P (perm)=0.001), except kelp whose composition was similar to that of any of the four other functional groups.

Figure 4 (A) PCO plot (based on Bray-Curtis similarity matrices) of the 12 fatty acids exhibiting at least 70% correlation in the six animal species, two macroalgal species, and two environmental components (see Table 1 for species list) sampled inside or outside of the South site (see Figure 1A). (B) Typical fatty-acid trophic biomarkers for those fatty acids included in the analysis (adapted from Parrish (2013) and Legeżyńska et al. (2014).

Of the three essential FAs (EPA, DHA [docosahexaenoic acid], and ARA), EPA was the most prevalent, present in all food web components except kelp and particularly abundant among the six animal species (Table S.1). ARA was in all components except seston, peaking in A. rubens (20%) and S. droebachiensis (15%) (Table S.1). Together, EPA and ARA contributed to 46% and 62% of the similarities in S. droebachiensis and A. rubens diets, respectively (Table B.1). DHA was present in all components except kelp, and was nevertheless less abundant than EPA and ARA, peaking at 9% in H. arctica (Table S.1).Animal and kelp FA profiles were generally dominated by polyunsaturated FA (PUFA), which ranged from 44% in Tonicella spp. to 65% in S. droebachiensis, and to a lesser extent by monounsaturated FA (MUFA), which varied from 17% in L. digitata and S. droebachiensis to 30% in Tonicella spp. (Table 3).

Animals and kelp contained lower levels of saturated FA (SFA), with lowest and highest proportions in respectively A. rubens (11%) and O. aculeata (27%) (Table 3). Conversely, seston and sediment contained more MUFA (37.9% and 42.5%, respectively) than PUFA (20.2% and 25.9%, respectively) (Table 3). SFA levels were higher in seston (41.6%) than sediment (24.2%) (Table 3). Animals exhibited the highest ratio of polyunsaturated to unsaturated FA (P/S; 3%), followed by kelp (2%), and seston and sediment (1%) (Table 3). All components, except for seawater, had a higher proportion of ω3 (omega-3) FAs than ω6 (omega-6). Animals and kelp contained about 3 times as many ω3 FAs (37.1% and 29.7%, respectively) than ω6 (11.9% and 10.5%, respectively) and sediment about twice as many ω3 FAs (13.0%) than ω6 (6.5%), while seston contained about half as many ω3 FAs (2.6%) than ω6 (5.3%) (Table 3). The DHA/EPA ratio was highest in seston (0.9), intermediate in sediment, A. rubens, H. arctica, and S. droebachiensis (0.1 to 0.5), low in Nereis spp., O. aculeata, and Tonicella spp. (0.02 to 0.03), and null (0) in kelp (Table 3).

Table 3 Sample size (N), mean proportional sum (Σ) of saturated (SFA), monounsaturated (MUFA), polyunsaturated (PUFA), ω3 (omega-3), and ω6 (omega-6) fatty acids, and mean ratios of polyunsaturated:saturated (P/S) and DHA EPA (DHA/EPA), in the six animal species, two macroalgal species, and two environmental components (see Table 1 for species list) sampled inside (I) or outside (O) of the South site (see Figure 1A).

3.4 Stable Isotopes and Trophic Magnification

Stable carbon isotope ratio (δ¹³C) differed significantly among the ten food web components included in the carbon isotope analysis (PERMANOVA, Pseudo-F_9,26 ⁼ 40.241, P (perm)<0.001), ranging from most depleted in seawater (-26.6‰) to least depleted in L. glaciale and A. rubens (-18.9‰) (Table 4). The δ¹⁵N values, which were lowest in L. digitata (3.4‰) and highest in A. rubens (11.0‰) (Table 4), also differed significantly among the nine food web components included (i.e. all components except Nereis spp.; PERMANOVA, Pseudo-F_8,25 ⁼ 130.64, P (perm)<0.001), indicating distinct trophic levels (see below). Hierarchical clustering analysis of δ¹³C and δ¹⁵N separated the latter nine components in four distinct groups (PERMANOVA, Psuedo-F_3,22 ⁼ 53.25, P (perm)<0.001; Figures 5 and S.1). Two of these groups each contained all samples of a single food web component, namely seston and L. digitata (kelp), hereafter termed respectively Group 1 and Group 4 (Figure 5). Group 3 contained three subgroups, each also containing all samples of a single food web component: sediment (Group 3a), L. glaciale (rhodolith, Group 3b), and H. arctica (Group 3c) (Figure 5). Group 2 had four subgroups, of which two were monospecific: A. rubens (Group 2a) and O. aculeata (Group 2d), and two each contained two species; A. rubens and Tonicella spp. (Group 2b) and S. droebachiensis and Tonicella spp. (Group 2c) (Figure 5). Group 1 (seston) had a significantly lower δ¹³C than all other groups (Tukey HSD, p<0.001), yet its δ¹⁵N was similar to that of Group 3 (infauna, L. glaciale, H. arctica) (Tukey HSD, p=0.779) (Figure 5). Group 2 (A. rubens, O. aculeata, Tonicella spp., and S. droebachiensis) had a significantly higher δ¹⁵N than Group 1, Group 3, and Group 4 (L. digitata) (Tukey HSD, p<0.001). Group 2’s δ¹³C was also significantly more enriched than that of Group 1 (Tukey HSD, p<0.001) and Group 3 (Tukey HSD, p=0.021), but not Group 4 (Tukey HSD, p<0.062). Group 4 had a significantly lower δ¹⁵N than all other groups (Tukey HSD, p<0.001) (Figure 5).

Table 4 Sample size (N), bulk stable isotope ratio (δ¹³C and δ¹⁵N; ‰), and relative trophic position (TP) in the six animal species, two macroalgal species, and two environmental components (see Table 1 for species list) sampled inside (I) or outside (O) of the South site (see Figure 1A).

Figure 5 Biplot of bulk carbon (δ¹³C) and nitrogen (δ¹⁵N) stable isotope ratios of five animal species, two macroalgal species, and two environmental components (see Table 1 for species list) sampled inside or outside of the South site (see Figure 1A). Nereis spp. was not included because of insufficient tissues for quantification in the N analysis. Components grouped (circled) based on agglomerative hierarchical cluster analysis (see Figure S.1).

At our assumed Δ¹⁵N fractionation factor of 3.4‰ (see section 2.8), the 10 food web components encompassed over three trophic positions (TP), with kelp (TP=1) and A. rubens (TP=3.2) at the base and top of the web, respectively (Table 4). Seston, sediment, rhodoliths (L. glaciale), and H. arctica occupied intermediate positions ranging from 1.8 to 2.1, whereas the three remaining animal species sampled had similarly high positions of 2.6 to 2.9 (Table 4). Thirty-seven (37) FAs correlated significantly with δ¹⁵N. Of those FAs, four (16:4ω3, 20:1ω11, 20:4ω6 [ARA], and 20:5ω3 [EPA]) exhibited a TMF > 1, and hence were biomagnified through trophic levels, whereas 33 had a TMF < 1, indicating biodilution (Table 5).

Table 5 Trophic multiplication factor (TMF) of 37 fatty acids (FA) as calculated from the slope (m) of corresponding linear relationship between FA concentration and bulk nitrogen (δ¹⁵N) stable isotope ratio (see section 2.8).

4 Discussion

Our study of the trophodynamics of a Newfoundland rhodolith bed is the first attempt to characterize nutritional patterns and trophic linkages of a rhodolith bed community with combined use of lipid, FA, and stable isotope analyses. Isotope-based rhodolith studies in the Northeast Atlantic (Grall et al., 2006) and California (Gabara, 2014) suggest macroalgae-based detritus are a key food source within rhodolith communities. Our comprehensive approach supported these findings, while showing that kelp may not be a significant food source among the 10 rhodolith bed components we considered. We identified three distinct trophic levels - producers, suspension/filter feeders and grazers, and predators - and discovered a potentially specific link between a macroalga-based diet and carbon source in mollusks and ophiuroids.

4.1 Rhodolith Community

Rhodolith communities vary globally in terms of numerically dominant cryptofauna, with gastropods dominating in the Maltese Islands (Sciberras et al., 2009), crustaceans in Santa Catalina Island (Gabara, 2014) and Ireland (de Grave and Whitaker, 1999), echinoderms and annelids in the Gulf of California (Foster et al., 2007), and polychaetes in South Australia (Harvey and Bird, 2008). In Newfoundland, chitons (Tonicella marmorea and T. rubra) and daisy brittle star (Ophiopholis aculeata) are the most common rhodolith epifauna (present study, Gagnon et al., 2012; Bélanger, 2020). The preponderance of chitons and daisy brittle stars at our study site, as well as overall species composition, aligned with the first detailed rhodolith bed biodiversity report from Newfoundland and Labrador (Gagnon et al., 2012). The consistency of chiton and brittle star abundances in conjunction with similar overall species composition in Newfoundland after over eight years, suggest high community stability within this bed. This notion is also supported by highly similar rhodolith morphological traits (shape and size) and biomass (19.5 kg m^-2) between both studies (present study, Gagnon et al., 2012). In the present study, most observed rhodolith cryptofauna was juvenile-sized, supporting the view that rhodolith beds are nursery grounds for a number of marine invertebrates (Foster, 2001; Kamenos et al., 2004b; Steller and Cáceres-Martínez, 2009; Gagnon et al., 2012).

4.2 Lipid Content and Classes

Lipid structure can vary based on environmental conditions, food availability, metabolism, and reproductive strategies (Fraser, 1989; Lloret and Planes, 2003; Parzanini et al., 2018). Low temperatures affect organisms by inducing changes in cell membrane fluidity, structure, and function (Crockett, 1998; Parrish, 2013; Colombo et al., 2017). To accommodate cold temperatures, cold-water ectotherms exhibit homeoviscious adaptation, a process of reducing sterol molecules and lengthening and unsaturating phospholipids (Hazel et al., 1991) to change lipid structure (Hall et al., 2000; Copeman and Parrish, 2003; Parrish, 2009). In the present study, we showed that the lipid structure of all samples varied by component, with overall high proportions of phospholipids (45 – 76%) and unsaturated FAs (58 – 87%), low proportions of sterols (6 – 19%), and a correlation between total lipids and sterol proportions. These findings strongly suggest increased cell membrane fluidity, while supporting our hypothesis that the lipid composition of our study organisms generally reflects the predominantly cold-ocean conditions of Newfoundland (0.3°C in April in the present study).

Triacylglycerols are a key component of lipid structure and the primary energy storage molecules. They are important to organisms during stressful periods such as limited food availability or reproduction. High variability of TAG content among organisms can result from differences in allocation strategies (i.e. reproduction, growth, or survival) (Fraser, 1989; Lee et al., 2006; Parzanini et al., 2018). As such, elevated TAG levels in H. arctica (Lebour, 1938), Tonicella spp., and O. aculeata (Himmelman et al., 2008) likely correspond to an abundance of food (i.e. phytoplankton) or reproductive timing; organisms increase and maintain their energy storage as they prepare for the reproductive season (Vanderploeg et al., 1992). Given the timing of our study (during the April phytoplankton bloom), high levels of TAG in H. arctica (34%), Tonicella spp. (29%), and O. aculeata (25%) likely demonstrate this energy storage trend and potential benthic-pelagic coupling, utilizing phytoplankton from the water column in the benthos as it becomes seasonally available (Iken et al., 2001; Hacker Teper, 2022). However, low TAG levels and high PL and ST levels of predatory/omnivorous A. rubens (2%, 76%, 10%, respectively), Nereis spp. (4%, 63%, and 19%), and S. droebachiensis (6%, 61%, and 15%) probably link to organisms in search of food for rapid growth during harsh winter conditions (Luis and Passos, 1995; Lee Jr et al., 2006). Over prolonged periods, low TAG levels could indicate stress (Fraser, 1989).

While seemingly low, lipid content levels of animal samples (~1% g g^-1 wet weight) were similar to those in other studies (Allen, 1968; Parzanini, 2018). In comparison to the animals, seston samples were mainly comprised of PL (50%) and FFA (31%), with no evidence of TAG (0%). This finding challenges previous reports of Newfoundland seston rich in AMPL and low in FFA (Parrish et al., 1995). Although rare, high levels of FFA could reflect degradation of lipids as a result of sewage runoff from St. Philip’s township (Parrish et al., 1992; Galois et al., 1996). Despite the unusual lipid profiles of our seston lipids, total lipid levels (~57 µg L^-1) were comparable to samples collected from a nearby Newfoundland site in 1991 (Parrish et al., 1995), likely evidence of the annual phytoplankton bloom with fresh, lipid rich diatoms (Budge and Parrish, 1998; Kiriakoulakis et al., 2005; Parrish et al., 2005). A sample volume greater than the 3 L of seawater filtered likely would have helped lipid class analyses. We therefore recommend filtering at least 10 L of seawater for future studies. Total lipid content of L. digitata samples (~1.3 mg g^-1 wet weight) were markedly lower than the wide range conveyed in the literature (~11 – 60 mg g^-1) (Raven et al., 2002; Schaal et al., 2010). A small portion of this difference could be because of the different techniques employed; gravimetric assays are typically 10 to 15% higher than Iatroscan-derived lipids (Parrish, 2013). Such low lipid contents could also reflect chemical changes of L. digitata in response to sea-ice induced environmental changes. However, brown seaweeds in subarctic cold waters have significantly more total lipid than those in tropical warm waters (Terasaki et al., 2009; Nomura et al., 2013), and brown seaweeds increase their total lipids in winter and under low light conditions (Honya et al., 1994; Nelson et al., 2002; Nomura et al., 2013). This pattern suggests our noted low total lipid content in L. digitata could indicate issues with the extraction of L. digitata lipid.

54.3 Fatty Acids and Stable Isotopes

Higher levels of unsaturated FAs compared to saturated FAs typically result from cold-water conditions (Parrish, 2009). However, unlike lipid structure, which organisms solely regulate (Arts et al., 2009), FA composition also depends on diet, feeding strategy, and phylogeny (Dalsgaard et al., 2003; Makhutova et al., 2011). High levels of ω3 FAs eicosapentaenoic acid (EPA; 20:5ω3) and docosahexaenoic acid (DHA; 22:6ω3), and occasionally ω6 FA arachidonic acid (ARA; 20:4ω6), typically characterize the marine environment. Our results suggest that except for H. arctica (9%) and A. rubens (5%), only trace amounts of DHA occur within the Newfoundland rhodolith community. Hiatella arctica and A. rubens likely require more DHA because mollusks need it for growth (Wacker et al., 2002; Arts et al., 2009) and hatching in copepods (Arendt et al., 2005), which were abundant in the diet of A. rubens. The relatively high levels of DHA in H. arctica mirror findings of Copeman and Parrish (2003) who reported that bivalves conserve relatively higher levels of plankton-derived DHA, and lower ARA, than echinoderms. The enhanced trophic relationship (TMF > 1) between A. rubens and calanoid copepod-derived lipids such as 20:1ω9 and 20:1ω11 (TMF, 1.19), whether feeding directly or indirectly, follows the suggestion from Connelly et al. (2014) that organisms can maintain lipid-rich energy stores from copepod-derived lipids in the same way as from DHA. The trophic relationship also aligns with our diatom abundance findings because calanoids in Newfoundland also rely on a diatom-based diet in winter (Urban et al., 1992; Beaugrand et al., 2002). Based on its nitrogen (δ15N) signature (11.0) and TP (3.2), we identified A. rubens as the top consumer of the studied components in the rhodolith food web. While we cannot ascertain the sources of its diet, its similar FA composition to H. arctica, combined with the predatory lifestyle of A. rubens towards mollusks (Allen, 1983) suggests H. arctica may be a potential prey item.

Eicosapentaenoic acid was the most abundant FA in each of the six animals sampled (A. rubens, 32%; H. arctica, 19%; Nereis spp., 29%; O. aculeata, 27%; S. droebachiensis, 25%; and Tonicella spp., 22%), which suggests a compensatory role for DHA deficiency. High levels of essential FA typically reflect the dominant microalgal group (Dalsgaard et al., 2003); EPA levels are consistently high in diatom-dominated environments, whereas DHA prevails where dinoflagellates dominate. We conducted our study at the beginning of the spring phytoplankton bloom, so high concentrations of EPA (diatoms) might have overshadowed any presence of DHA, such as in H. arctica and A. rubens (Budge et al., 2001; Dalsgaard et al., 2003). As indicated by its prevalence in all diets, EPA has a high trophic magnification factor (TMF, 1.15). The biomagnification of EPA may confirm our hypothesis of a bottom-up food web in which most organisms rely on a shared resource (diatoms, EPA) passing from first order consumers onto second and third order consumers. The presence of diatom trophic biomarkers 16:0 and 16:1ω7 in sediment and in all animal diets, further points to diatoms as a major food source in this rhodolith community. Although synthesized biomarkers like 16:0 are less useful for understanding dietary intake than externally derived FAs like EPA, we recommend considering synthesized FA given their obvious accumulation in the present study (Wennberg et al., 2009). Interestingly, seston samples contained little EPA (1.2%) or 16:1ω7 (1.2%), but rather large proportions of 16:0 (21.7%; diatom and particulate macroalgae), 18:1ω9 (32.9%; zooplankton and particulate macroalgae), and 18:0 (18.2%; detrital) biomarkers (Wakeham and Canuel, 1988; Parrish, 2013; Legeżyńska, et al., 2014). This is in accordance with Bec et al. (2010) who concluded phytoplankton only explains 27% of the variance in seston. Because of our study site’s proximity to riverine input, some proportions of 18:1ω9 and 18:2ω6 (5.3%) could also be influenced by conifer pollen found in large quantities during spring (Masclaux et al., 2013; Lichti et al., 2017). Regardless, diatoms help support a rich infaunal community living underneath rhodoliths, which likely depends on deposition of organic material at the rhodolith-sediment interface (Steller et al., 2003; Grall et al., 2006; Berlandi et al., 2012). However, diatoms may provide more than just food to rhodolith beds, as per Steller and Cáceres-Martínez’s (2009) suggestion that diatom films on rhodoliths promote larval settlement of invertebrates (Morse et al., 1988; Daume et al., 1999; Huggett et al., 2006).

Although less abundant than EPA or DHA, at least some marine organisms also require ARA as an essential FA. ARA is particularly important for echinoderms (Copeman and Parrish, 2003) to regulate metabolic activities (Ciapa et al., 1995) and to maintain membrane structure and function (Parrish, 2009). In our study, ARA was highest in the echinoderms A. rubens (20%) and S. droebachiensis (15%). These high levels of ARA in sampled urchins were higher than those reported by Kelly et al. (2008) from both a coralline barren and kelp bed (Laminaria digitata). Laminaria digitata, a common food source for S. droebachiensis in eastern Canada, is typically rich in ARA (Fleurence et al., 1994; Schmid et al., 2014) and its low nitrogen (δ¹⁵N) signature (3.4) places it at the base of our food web (TP: 1.0). However, our kelp samples contained only 0.1% of ARA proportionally, making it an unlikely source of ARA across our food web. Nevertheless, high levels of a FA precursor to ARA, 18:2ω6 (9%), marked FA profiles of L. digitata. The ability of S. droebachiensis to synthesize ARA from 18:2ω6 complicates biomarker identification from algal FAs (Kelly et al., 2008). Combining FA analyses with stable isotope analyses can increase the resolution of organism diets. Thus, given some similarities of carbon (δ¹³C) signatures, S. droebachiensis (-19.8‰) may consume L. digitata (-21.0‰); however, δ¹³C signatures of L. digitata were degraded compared to literature (Raven et al., 2002; Schaal et al., 2010). Although we could not identify the dietary source of ARA among sampled components, we conclude amphipods represent a potential source of ARA in the rhodolith system because of both their high abundance in the rhodolith community and their typical ARA richness (Guerra-García et al., 2004). Given its biomagnification across diets (TMF, 1.16), ARA likely represents a key essential FA in the rhodolith food web and should be a focus of future studies.

Three potential photosynthetic carbon sources exist within a rhodolith bed: macroalgae (including the rhodoliths themselves), phytoplankton, and microphytobenthos (Grall et al., 2006). Similarly to S. droebachiensis, the (δ₁₃C) signature (-19.3‰) of Tonicella spp. resembled that of rhodoliths (-18.9‰) and not kelp, suggesting that our focal grazing species feeds directly on coralline algae with little to no dietary input from kelp. Strongylocentrotus droebachiensis often leaves star-shaped tooth marks on rhodoliths (personal observations) and consumes coralline algae and microalgal films in coralline barrens, though limiting their bites to tips of rhodolith branches (Steneck, 1990; Scheibling et al., 1999; Lauzon-Guay and Scheibling, 2007). Tonicella spp. occurs in high abundance, often on the outside of rhodoliths. Although their bite marks are not as evident as those of S. droebachiensis, their articulating plates give them the unique ability to graze between rhodolith branches (Steneck, 1990). In addition to consuming corallines, these grazers can graze on diatom films on the surface of rhodoliths (James, 2000). Lipid concentrations were too low to confidently rely on our rhodolith lipid and FA results, therefore preventing comparisons with results for S. droebachiensis and Tonicella. spp. Contrary to Kelly et al. (2008) who sampled a blend of coralline algae taxa, our study is the first to explore lipids in a (presumably) single rhodolith species (L. glaciale). We were uncertain of the quantity of rhodolith material needed for our lipid analyses. We found that 3-5 g of rhodolith material yields negligible lipid concentrations. Based on our subsequent calculations and recent findings by Teed et al. (2020) that Newfoundland rhodoliths are 85% CaCO₃ by weight (Teed et al., 2020), we recommend a minimum of 12-15 g of rhodolith material (organic and inorganic, combined) for future analyses. As mentioned above, δ¹³C signatures of L. digitata were lower than published accounts. We propose that sea ice contributed to carbon depletion in this species by causing bacterial degradation of kelp tissues or the production of new fronds on stored photosynthates containing little nutritional value, or by creating inadequate phytosynthetic conditions (Fredriksen, 2003; Vanderklift and Bearham, 2014). Furthermore, because δ₁₅N signatures are influenced by light availability (Vanderklift and Bearham, 2014) and nutrient status (Gagné et al., 1982; Schaal et al., 2009), it is possible that sea ice also contributed to the unusually low δ¹⁵N signatures of L. digitata (3.4‰) (Raven et al., 2002; Schaal et al., 2010). This and the low total lipids in L. digitata require additional investigation.

Fatty acid profiles of H. arctica, O. aculeata, and Tonicella spp. were influenced by bacterial (16:2ω4 and i17:0) and red and green seaweed (18:1ω7) trophic biomarkers. δ¹⁵N signatures identified H. arctica (6.3‰) as a first order consumer and O. aculeata (9.0‰) and Tonicella spp. (9.7‰) as second order consumers. Though their trophic levels differed (TP, 1.9, 2.6, and 2.9, respectively), the similar δ¹³C signatures (-19.3 – -20.0‰) in the three species indicate similar carbon sources. In combining FA and stable isotope analyses, we identified similarities between the diets of non-predatory benthic organisms that suggest a resource partitioning relationship wherein animals consume different elements of the same foods (i.e. particulate algae, microphytobenthos, phytoplankton) at different times using different strategies (filter feeding, suspension feeding, and grazing, respectively), resulting in the effective use of the majority of the food source (Hines, 1982; Parrish et al., 2009). Accordingly, the feeding strategy of O. aculeata offers a possible explanation for the feeding relationship described above. Although brittle stars shared a nearly identical δ¹³C signature to H. arctica (-19.9‰ and -20.0‰, respectively), their stronger δ¹⁵N signature (9.0‰ versus 6.3‰, respectively) presumably originated from POM enrichment as POM settled on the seafloor (Iken et al., 2001). Therefore, O. aculeata likely consumes the same material as H. arctica, but as resuspended matter rather than through direct filter feeding, and thus benefits from the particulate leftover from H. arctica. Additionally, Tonicella spp. may graze on larger particulate material inaccessible to O. aculeata, potentially breaking it down into smaller pieces for consumption by O. aculeata. Indeed O. aculeata utilizes several different suspension feeding mechanisms to collect from the water column enriched, resuspended benthic particles (Labarbera, 1978). The species may therefore enhance benthic-pelagic coupling in transferring organic material from the pelagic zone to the benthos (Guerra-García et al., 2004). Our findings corroborate Grall et al.’s (2006) demonstration that rhodolith beds comprise both pelagic (Herman et al., 2000) and benthic, algae-based feeding relationships (Takai et al., 2004).

4.4 Conclusion and Future Research Directions

Our study supports the hypotheses that (1) the lipid composition of rhodolith cryptofauna generally reflects the predominantly cold-ocean conditions of Newfoundland; and (2) the rhodolith bed food web is mainly controlled from the bottom up by planktivores and detritivores as reflected by high abundance of planktonic and bacterial biomarkers. Our lipid and FA analyses revealed high levels of phospholipids and unsaturated FAs combined with low sterols in all animal species, suggesting adaptability for enhanced cell membrane fluidity in a cold-water environment. Our FA and stable isotope analyses showed that many taxa in the rhodolith community rely on a shared resource: diatoms. We also unveiled potential resource partitioning between first- (H. arctica) and second- (O. aculeata and Tonicella spp.) order consumers, whereby differences in feeding strategies enable utilization of specific components of the same organic and inorganic material.

Our study documents, for the first time, the specific diets of, and feeding relationships among, dominant animal taxa in a cold-ocean, Newfoundland rhodolith bed. The complexity of feeding relationships in benthic systems (Kharlamenko et al., 2001; Pitt et al., 2009; Kelly and Scheibling, 2012), temporal variability of benthic-pelagic coupling (Iken et al., 2001; Hacker Teper, 2022), and uncertainty in the amount of sample material required for some bed components (seawater and L. glaciale), limited the ability to decipher all possible linkages among and above the ten rhodolith bed components included in the present study. Future studies of rhodolith bed trophydynamics should use the aformentioned minimal sampling material requirements. We also recommend expanding the breadth of focal species to include, for example, pelagic organisms such as amphipods and copepods (Guerra-García et al., 2004; Pakhomov et al., 2004). The unknown contribution to trophodynamics of sediment and associated infauna underlying rhodoliths should be explored, particularly that of bacteria and microbes and their production of FAs (e.g. 16:0, 16:1ω7, and 18:1ω7) (Fullarton et al., 1995). The present study provides a snapshot of the feeding relationships in a subartcic rhodolith bed. A companion study at our study site suggests the feeding relationships reported in the present study are spatially and temporally stable with no appreciable intra annual changes outside of the spring phytoplankton bloom (Hacker Teper, 2022). Similar studies with broader spatial and temporal coverage would help further understand the diversity and stability of these relationships in rhodolith communities.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material. Inquiries can be directed to the corresponding author.

Author Contributions

SHT: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Field Work, Data curation, Writing - original draft, Visualization, Project administration. CP: Conceptualization, Methodology, Investigation, Resources, Writing - original draft, Supervision, Project administration, Funding acquisition. PG: Conceptualization, Methodology, Investigation, Field work, Resources, Writing - original draft, Supervision, Project administration, Funding acquisition.

Funding

This research was funded by the Natural Sciences and Engineering Research Council of Canada (NSERC Discovery Grants to C. Parrish: RGPIN-2017-04639 and 16940: RGPIN-2015-05343), the Canada Foundation for Innovation (Leaders Opportunity Fund Grant to 16940), and the Research & Development Corporation of Newfoundland and Labrador (Ignite R&D Grant to 5003.070.002).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We are grateful to Jeanette Wells and Alison Pye for assistance with laboratory work, David Bélanger, Anne Provencher St-Pierre, Samantha Trueman, and the Ocean Sciences Centre Field Services Unit for assistance with field work, and Paul Snelgrove for providing access to CTD equipment. We also thank Paul Snelgrove, Suzanne Dufour, and two anonymous reviewers for constructive comments that helped improve the manuscript.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmars.2022.899812/full#supplementary-material

References