In addition to its extreme dryness, the Atacama Desert holds records for other extreme environmental characteristics. The highest surface ultraviolet (UV) radiations and total solar irradiances ever measured on Earth have been reported from its high altitude "Altiplano" area (Piacentini et al., 2003; Cordero et al., 2014) and pre-Andean Domeyko Cordillera (Rondanelli et al., 2015). The prevalent cloudless conditions and relatively low total ozone column in the Atacama Desert result in extremely high solar radiation, increasing considerably the probability of lethal photoinhibition and photooxidative damage to phototrophs (Solovchenko and Merzlyak, 2008). Adaptation strategies for protection against high solar irradiance were observed in several microorganisms from the Atacama Desert. For example, the thallus of the lichen Acarospora schleicheri, which is found in the Altiplano area of Chile, accumulates rhizocarpic acid upon exposure to high levels of UV-B radiation (Rubio et al., 2002). Biomineralization of amorphous silica (sinter) around cyanobacteria filaments has been reported as a significant irradiation shield from the El Tatio geothermal field at 4300 m a.s.l. in the Andes Mountains (Phoenix et al., 2006). A different strategy to avoid stressful light conditions was found within the endolithic habitat of halite (NaCl) rocks.

When both extreme desiccation and high ambient UV fluxes are combined, such as in areas of the hyperarid core of the Atacama Desert, the surface of rocks can be rendered sterile (Cockell et al., 2008). Under these extreme conditions, a millimeter-thick layer of gypsum has been found sufficient to protect microorganisms from UV-induced killing (Cockell et al., 2008).  In many of the above studies, translucent gypsum was reported to allow photosynthesis while providing protection from excess solar irradiance and retaining sufficient moisture for microorganisms survival. These factors may account for the high biodiversity found in gypsum-hosted communities (see Supplementary Table S1 for a summary of gypsum endolithic communities). However, there is a lack of knowledge regarding the microorganisms vertical organization or the distribution of photosynthetic pigments and secondary metabolites within endolithic community. No studies have addressed the effect of light gradients, or other characteristics of the lithic substrate, on the spatial localization of the different microbial taxa inhabiting this extreme environment.

The aim of the present work was to elucidate exceptional adaptation strategies of chlorophototrophic (i.e., chlorophyllbased phototrophs) microorganisms colonizing the insides of gypsum deposits in the hyperarid zone of the Atacama Desert. These endolithic microbial colonies are taking advantage of the physical and structural properties of translucent gypsum deposits containing sepiolite inclusions. Using a combination of microscopic and spectroscopic analytical methods, molecular and phylogenetic analyses, and micromorphological studies, we provide a detailed characterization of the chlorophototrophic microorganisms colonizing the gypsum endolithic habitat. Together with long-term microclimatic data and the composition and structure of the mineral substrate, we constructed a conceptual model of the vertical succession of microorganisms, together with the distribution of their photosynthetic pigments and secondary metabolites, within gypsum deposits under extremes of aridity and solar radiation.

# Materials and Methods

The following is the summary of methods used in this study. More detailed information is provided in Supplementary Materials.

# Site Characterization and Sampling

The sampling zone (23◦ 53′ S, 068◦ 08′W and 2720 m a.s.l.) was located in the south of the Salar de Atacama basin (**Figure 1A**) in colonization.

the Socompa volcano is seen on the horizon. (C) Cross section of gypsum bearing cryptoendolithic (red and green arrows) and hypoendolithic (blue arrow) microbial

the north-south trending depression of the Cordon de Lila range, in northern Chile. This depression is mostly covered by volcanic ash material but large gypsum outcrops (1 × 2 m) can be found in several locations (**Figure 1B**). In the field, gypsum deposits were assessed for endolithic colonization by visual inspection of microbial pigments present in fractured samples. These pigments indicate the presence of cryptoendolithic (occupying pore spaces beneath rock surface, Nienow, 2009) and hypoendolithic (colonizing the undermost layer of the rock, sensu Wierzchos et al., 2011) microbial colonization forming horizons beneath the surface and close to the bottom of the gypsum deposits (**Figure 1C**). Colonized gypsum samples (10 × 10 cm) were collected in 2013 in the Cordon de Lila area, sealed in sterile Whirlpacks R , and stored at room temperature in a dark and dry environment until analysis.

# Environmental Data Acquisition

Microclimate data were collected 6 km from the sampling site for 39 months (37 for T and RH), from January 20th, 2010 to April 20th, 2013, using an Onset HOBO R Weather Station Data Logger (H21-001) connected to a SolarStream R solarpowered transmitter for data transmission by Iridium Satellite Constellation. Sensors recorded measurements every 30 min intervals. Air relative humidity (RH) and temperature (T) were recorded 25 cm above the rock surface (probe was shaded from the sun) using RH/T sensors (HOBO R S-THB-M002; precision, ±2.5% RH and ±0.2◦C T). Solar flux was measured using a photosynthetically active radiation (PAR) sensor for wavelengths of 400–700 nm. PAR data were also indicative of cloud cover at the sampling site. Rainfall was monitored using a Rain-o-Matic 100 (PRONAMIC ApS, Herning, Denmark) tipping bucket gauge (resolution of 1 mm). The presence of liquid water on the gypsum surface was determined by electrical conductivity (EC) using sensor as described by Wierzchos et al. (2013) (see also Supplementary Section S1).

# Mineralogy X-ray Diffraction Analyses

The mineralogical composition of the gypsum duricrust was identified by X-ray powder diffraction (XRD) using a Philips X'Pert diffractometer (Philips, Amsterdam, The Netherlands) with graphite-monochromated CuKα radiation (see detailed procedures in Supplementary Section S2).

# Petrographic Microscopy

Petrography studies of gypsum duricrust thin sections (30µmthick) were conducted using a Nikon Eclipse LV100Pol (Nikon, Tokyo, Japan) polarized light microscope equipped with a Nikon DS-Fi1 digital camera.

# Mercury Intrusion Porosimetry

Pore space and pore size distribution of the gypsum samples were characterized by mercury intrusion porosimetry (MIP). Three types of subsamples were analyzed: samples from the surface dense crust (above the colonization zone), samples from the interior of the gypsum deposits (within the colonization zone), and samples from the bottom dense crust (covering the hypoendolithic colonization zone). A PoreMaster 60/Quantachrome (Quantachrome Instruments, Boynton Beach, USA) instrument was used to determine porosity (connected porosity) of the rock samples in the pore diameter range 0.0036–190.58µm. A mercury surface tension value of 0.48 N/m and rock-mercury contact angle of 140 dag were used in the Washburn or Laplace equations, respectively.

# Measurements of Light Conditions within Gypsum Deposits

The photosynthetic photon fluence rate (PPFR) was measured beneath the gypsum surface at different depths (up to 1.5 mm) with a small (3.5 mm in diameter) spherical sensor. These measurements were performed on wetted gypsum after simulation of a 10 mm rainfall to measure light conditions, inside the rock, when the presence of water would allow photosynthetic activity. Measurements for UVA and UVB radiation beneath the gypsum surface at different depths (up to 1.5 mm) were performed using cosine corrected sensors. Detailed procedures are provided in Supplementary Section S3.

# Fluorescence Microscopy (FM)

Small fragments of gypsum showing distinct signs of endolithic colonization with a green/orange layer beneath the rock surface and a green layer just above the bottom of the gypsum sample were scraped and suspended in double-distilled water. Suspensions were vortexed (2 min) and allowed to settle for 1 min before removing the supernatant, which was centrifuged for 10 min at 5000 g. The pellet was then stained with SYBR Green I (SBI) (Molecular Probes), a fluorocromes used for the specific staining of nucleic acids. Bright field images, SBI fluorescence and photosynthetic pigment autofluorescence were visualized in a Zeiss AxioImager D1 fluorescence microscope (Carl Zeiss, Jena, Germany). Detailed procedures are provided in Supplementary Section S4.

# Scanning Electron Microscopy in Backscattered Electron Mode (SEM-BSE)

Colonized gypsum samples were processed for SEM-BSE and energy dispersive X-ray spectroscopy (EDS) microanalysis as described in Wierzchos et al. (2011). As the intensity of the BSE signal depends on the mean atomic number of the sample, the SEM-BSE technique identifies heavy metal-stained ultrastructural elements of microbial cells. SEM-BSE was then used in combination with EDS to characterize the minerals associated with specific cell aggregates. Gypsum samples were observed using a scanning electron microscope (DSM960 Zeiss, Oberkochen, Germany) equipped with a solid-state, four diodes BSE detector plus an auxiliary X-ray EDS microanalytical system (Link ISIS Oxford, UK).

# Transmission Electron Microscopy (TEM)

Fragments of colonized gypsum material were processed for transmission electron microscopy. Ultrathin sections were visualized using EM910 (Leo, Oberkochen, Germany), CM200 Philips (Philips, Amsterdam, The Nederlands) and analyzed using energy dispersive X-ray spectroscopy (EDS) microanalysis. Detailed procedures are provided in Supplementary Section S5.

# Low Temperature Scanning Electron Microscopy

We characterize the features of gypsum deposits and internal micromorphology of microbial cells in their natural hydrated state using low temperature scanning electron microscope (LT-SEM) according to procedures described by Wierzchos et al. (2012b). The LT-SEM (DSM960 Zeiss, Oberkochen, Germany) was equipped with a cryotransfer system (CT1500, Oxford, UK) and with imaging systems using secondary electrons (SE) and backscattered electrons (BSE) detectors, and with an X-ray energy dispersive spectroscopy (EDS) system.

# Algal Culturing and Micromorphological Characterization

Fragments of gypsum samples bearing the green/orange colored cryptoendolithic colonization were visualized by stereoscopic microscope with aim to find algal colonies. Small rock fragments inhabited by algae were mechanically scrapped off into Petri dishes poured with agarized BBM (Deason and Bold, 1960) and BG-11 mineral media (Stanier et al., 1971). Algae were inoculated at 20◦C under a 14 h/10 h light/dark regime using photon irradiance of approximately 30–50µmol photons m−2 s −1 provided by 18 W cool fluorescent tubes (Philips TLD 18W/33). After 5–10 weeks, the agar plates were checked for the presence of algal microcolonies. The algae were investigated using an Olympus BX51 (Olympus, Tokyo, Japan) light microscope with differential interference contrast (DIC). Microphotographs were taken with an Olympus Z5060 digital camera.

# Microbial Community Analyses

Environmental DNA was extracted from the top colonization zone of two individual gypsum rocks using the PowerSoil DNA isolation kit (MoBio laboratories Inc., Solana Beach, CA) as described in Robinson et al. (2015). DNA was amplified using the barcoded universal primers 27F and 506R for the V1–V3 hypervariable region of the 16S rRNA gene and amplicons from 3 reactions were pooled together for sequencing using the Illumina MiSeq platform by the Genomics Resource Center (GRC) at the Institute for Genome Sciences (IGS), University of Maryland School of Medicine (Robinson et al., 2015). We obtained a total of 3326 paired-end reads for two rock samples (CL1 and CL2) collected from the same location, with an average size of 436 bp. The QIIME package (v1.6.0) was used to process the Illumina paired-end reads as previously described (Robinson et al., 2015). Sequence reads were normalized to 1662 reads per sample and diversity metrics were calculated based on OTUs at the 0.03% cutoff in QIIME (Caporaso et al., 2010). Nonparametric Mann-Whitney tests were used for statistical analysis. We use the following 18S rRNA primer sets to obtain nuclear sequences from the algae: 82F-1498R (Bachy et al., 2011); ChloroF-ChloroR and BaciF-BaciR (Valiente Moro et al., 2009); Euk328F-Euk329R (Romari and Vaulot, 2004); EUK528f-CHLO02r (Zhu et al., 2005).

# Analyses of Photosynthetic Pigments by Raman Spectroscopy

Cryptoendolithic and hypoendolithic colonization zones, harboring phototrophic eukaryotes and prokaryotes, were examined using point Raman analysis on a Thermo Fisher DXR Raman microspectrometer (Thermo Scientific, Waltham, USA) with a 532 nm diode laser and a Renishaw InVia Reflex Raman microspectrometer using both 785 nm diode and 514.5 nm Ar laser. The later Raman system was used to construct Raman distribution maps of carotenoids within the algal zone using 514.5 nm laser line. The analyses were performed on freshly cut transects of the gypsum substrate. The Raman spectroscopy technique was chosen for its in situ approach with small amount of biological material. Detailed procedures are provided in Supplementary Section S6.

# Spectrophotometric Analyses

The absorption spectra of acetone extracted photosynthetic and light-protecting pigments (chlorophyll a and b, and carotenoids) were obtained using a HP 8452A Diode Array dual-beam spectrophotometer (Hewlett-Packard, Tokyo, Japan). Absorbance values for characteristic wavelength (maximum peaks for particular pigments) were selected for semiquantitative determination of pigments contents [mg/g of the dry weight (DW) in powdered sample] using the trichromatic equations and extinction coefficients (Lichtenthaler, 1987). Detailed methodology sample processing and pigment extraction are provided in Supplementary Section S7.

# Results

# Microclimate Data

Microclimatic parameters were recoded with a microweather station in close proximity to the sampling site over 39 months (37 for T and RH measurements), starting in January 2010. The maximum daily photosynthetically active radiation, rainfall precipitations, and electrical conductivity measured on the gypsum gypcrete surface, as well as mean the daily air relative humidity and mean daily air temperature oscillation measured at the sampling site are shown in **Figure 2**. Some of this data (first 22 months) was previously reported by Wierzchos et al. (2013). PAR values revealed a very intense solar irradiance over the 3-year period of recording (**Figure 2A**). The mean daily PAR value was 1175µmol photons m−2 s −1 and during the austral summer (January, n = 4) the mean maximum daily value was as high as 2298µmol photons m−2 s −1 . For a number of days, usually in February, spikes of extremely high PAR values, up to 2554µmol photons m−2 s −1 , were recoded (**Figure 2A**). These spikes likely originated from light scattering by cumulus cloud formation (Piacentini et al., 2003), although cloudy days were very scarce.

The sum of recorded rainfall precipitation over the threeyear period was 88 mm and represented 19 separate episodes (**Figure 2A**). The annual mean precipitation was 27.1 mm y−1 with 61% of rain precipitations occurring during the summer months (January-March). The maximum daily rain precipitation was observed on 17th June 2011 with a total of 14 mm of rain between 8:30 a.m. and 11:00 p.m. EC sensors detected liquid water on the gypsum surface mostly in the late evening and until about midnight (**Figure 2A**). Most of these events corresponded to rainfall episodes with the exception of 7 events where moisture

was detected by EC readings with no rainfall gauge reading. These 7 events could be the result of intermittent rainfalls below the detection limit of rainfall gauge (<1 mm). Moisture on the rock surface persisted several hours after the rainfall (see Supplementary Figure S1). The total time of moisture at the rock surface, as detected by the EC sensor, was 638 h for the 39 months periods, giving a value of 196 h y−1 (about 8 days y−1 ). These values represent surface water and it is likely that moisture was available inside the gypsum rock for longer periods. Time periods between moisture events were exceptionally long with a maximum of 254 days in 2010. Analyses of our microclimate data indicated that rainfall was the only source of liquid water potentially available for the gypsum microbial communities.

Oscillation of mean daily air temperatures was in accordance with seasonal changes of solar irradiance (**Figure 2B**). The mean T for the 37 month period was 15.3◦C. The mean daily maximum T in summer (January-February 2013) was 43.6◦C with a maximum of 49.3◦C. The minimum recorded T was −7.4 ◦C. Air RH values were remarkably low over the entire recording period. Mean annual RH for 3 years was as low as 16.5% with mean winter and spring seasonal RH values between 7 and 11%. RH above 80–90% was only observed during short and infrequent rainfall episodes.

# Mineralogical and Physical Characterization of Gypsum-composed Gypcrete

Gypsum formations at the Cordon de Lila sampling site appear as hard layer deposits on the soil surface, interbedded between layers of ignimbrites and/or filling large cracks on these rocks. Volcanic ashes cemented by gypsum have been previously reported in this region of the Atacama Desert (Bao, 2005). This type of formation is called gypcrete according to the nomenclature by Horta (1980). The XRD analysis of Ca-sulfate deposits revealed gypsum (CaSO4·2H2O) as the major mineral (96–98%) with minor amounts of cristobalite (SiO2), calcite (CaCO3), and potassium feldspar (KAlSi3O8). A clay mineral with X-ray reflection at 12.2 Å was identified as sepiolite by in situ qualitative and quantitative EDS coupled to high resolution TEM (see Supplementary Figure S2). Study of thin sections with petrographic microscopy (results not shown) revealed dense morphologies with highly compacted large lenticular gypsum crystals in the surface layer (up to 2 mm deep) and in the bottom layer (2–3 mm) of the gypsum deposit. The appearance of these compact layers was indicative of gypsum dissolution and recrystallization processes. The interior of the sample consisted of microcrystalline (50µm long) lenticular gypsum crystals crossed by dense horizontal bands (0.3–1 mm-wide) of verticalized sub-prismatic gypsum crystals and lenticular gypsum crystals. Sequences of vertically-elongated pores (<2 mm-wide) and gypsum "columns" (0.1–1 mm-wide) were also observed within the interior of the gypsum deposits. These columns were formed by anhedral gypsum crystals and presented evident dissolution features. These petrological observations were consistent with porosity analyses performed with the MIP technique. MIP data showed a total porosity for connected pores (for a 0.0036µm to 190.58µm range) of 29.9% (n = 4, SD = 3.9) for the surface layer, 56% (SD = 10.1) for the rock interior, and 25.6% (SD = 3.8) for the bottom layer of the gypsum. Pore size distribution curves revealed larger pore diameters, between 0.1 and 190.58µm, for the gypsum interior than for the surface and bottom layers (**Figures 3A–A**′′). The total intrusion volume of mercury was about 0.08 ml g−1 for the surface and bottom layers and 0.31 ml g −1 for interior of the gypsum deposit. Together, these data support a significantly higher porosity for the interior than for the external layers of the gypcrete.

The transmitted PAR value within the gypsum colonization zone at a depth of 2–5 mm, we measured was 0.1–1% of that of the incident PAR light (**Figure 3B**). Using 2700µmol photons m−2 s −1 as the incident PAR value for a sunny day we calculated an absolute PPFR value between 0.2 and 30µmol photons m−2 s −1 within the cryptoendolithic colonization zone of the gypsum deposit. No UVA or UVB radiation was detected across a 1.5 mmthick gypsum layer.

# Characterization of Endolithic Microbial Communities

### Characterization of Endolithic Microbial Communities by Microscopy

In field fracturing of gypsum deposits revealed an abundant microbial colonization, which was clearly visible as a thick 2–6 mm orange-green layer beneath the gypcrete surface (**Figure 4A**). Closer examination of this cryptoendolithic habitat shows the presence of orange-colored algae cells close to the surface and green algae cells just below the orange layer of cells. This transition from green to orange cells was gradual as shown in **Figure 4A**′ . DIC microscopy of orange mature algae showed a cytoplasm filled with orange-colored lipids and small loose chloroplasts surrounding the nucleus (**Figures 4B**′ **,C,C**′ ). Chloroplasts within the green algae cells were denser and

deposits. (A–A**′′**) pore size distribution curves for surface layer (A), interior (A**′** ) and bottom layer (A**′′**); mean curves for four samples. (B) Intensity of photosynthetic photon fluence rate (PPFR) within gypsum wet deposits at different depths (n = 3); hypothetical position of cryptoendolithic colonization zone is shown by the schematic positioning of green-orange algae (circles) and cyanobacteria (small cyan discs).

found throughout the cytoplasm (**Figures 4D,D**′ **,E**). Green algae aggregate showed some cells with a lower density of chloroplast and starting to take orange color (**Figures 4E,E**′ ). In some places, it was also possible to observe cyanobacteria aggregates just beneath the algae layer (**Figure 4F**). Cyanobacteria aggregates were also found within the colonization zone at the bottom of the gypsum deposits (**Figure 4G**). This hypoendolithic habitat displayed a dense bottom layer and micro cave-shaped pores within the colonization zone (**Figure 4H**). Fluorescence microscopy revealed potentially viable cyanobacteria cells and cells with degraded photosynthetic pigments (Roldán et al., 2014) within this hypoendolithic habitat (**Figure 4G**).

The cryptoendolithic habitat appears beneath the 0.5–1 mmthin hardened surface layer with dense large lenticular gypsum crystals, as shown by the SEM-BSE technique (**Figure 5A**). In many places, large aggregates (0.5–1 mm) of clay mineral identified as sepiolite (SEM-BSE-EDS results not shown) were also observed (dash outline in **Figure 5A**). Algal cells were

FIGURE 4 | Cross-sections of gypsum layered microbial endolithic colonization. (A) Stereoscopic microscopy view with details in (A**′** ) of the orange-to-green cryptoendolithic colonization layer close to the gypsum dense surface (Surf); red, green, and blue arrows point to orange algae, green algae, and cyanobacteria colonization zones respectively. (B) DIC image of orange mature algae and (B**′** ) DIC image insert where some relicts of loose chloroplasts within the mature alga can be distinguish. (C) DIC and (C**′** ) FM image of orange algae with dispersed chloroplasts [red signal in (C**′** )] surrounding the alga nucleus [arrow in (C**′** )]. (D,E) DIC images and (D**′** ) FM image of green algal cells with uniformly distributed dense chloroplasts in the whole cytoplasm [red signal in (D**′** ); arrow in (D**′** ) point to nucleus and arrow in (E) point to chloroplasts; (E**′** ) DIC image of green cells in the first step of the carotenogenesis process with less dense chloroplasts (arrows). (F) MF image of cryptoendolithic cyanobacteria with red autofluorescence signal. (G) FM image of hypoendolithic cyanobacteria; autofluorescence red signal that potentially belong to vital cells and green signal corresponding to cells with degraded photosynthetic pigments. (H) Stereoscopic microscopy view of hypoendolithic cyanobacteria colonization zone (blue arrow) within gypsum close to the bottom (arrow) of the rock. Note the dense structure of the bottom (Btm) layer and small cave-shaped pores in the hypoendolithic habitat (asterisk). Scale bars = 20 µm, except (A,A**′** ,H).

observed between smaller gypsum crystals, just beneath the hardened surface layer. A detailed view of these algae with a cytoplasm filled by lipoidal globules is shown in **Figure 5B**. Algal cells found deeper in the cryptoendolithic zone contained fewer lipoidal globules and their cytoplasm was mostly filled by chloroplast structures (**Figures 5C–F**). Cyanobacteria cells were mostly found close to sepiolite aggregates (**Figure 5G**) beneath the green algae layer.

Similarly, to DIC microscopy analyses (**Figures 4B**, **6A**), visualization of algae from the orange zone with TEM showed abundant lipoidal globules occupying almost the entire cell (**Figure 6B**).  TEM of algae subjacent zone revealed cyanobacteria cells surrounded by electron dense structures, forming concentric sheaths of polysaccharides and abundant extracellular material embedding heterotrophic bacteria (**Figure 6C**). Our TEM

microscopic observations confirmed that the micromorphology of the gypsum cyanobacteria, in the top layer of the substrate, is similar to that of members of the Chroococcidiopsis genus (**Figure 6C**). Note that cyanobacteria from the cryptoendolithic habitat were always observed in close contact with floccules of sepiolite (**Figure 6C**). The mineralogical nature of this clay mineral was identified by in situ qualitative and quantitative EDS coupled to high resolution TEM (Supplementary Figure S2). **Figure 6D** shows Chroococcales cyanobacteria from the hypoendolithic habitat.

Visualization of cryo-fixed and cryo-fractured samples was essential in observing the gypsum microbial communities in situ in their "natural" state after sample rehydration (Ascaso and Wierzchos, 2002). It allowed for the characterization of internal structures in algae cells and that of the extracellular polysaccharidic substances (EPSs) surrounding microbial cells. LT-SEM of fractured algal cells from the green zone showed well defined chloroplasts filling almost the entire cytoplasm (**Figure 7A**). In contrast, algal cells from the orange zone accumulated lipoidal globules (**Figure 7B**). In many cases, we observed rod shaped (1–2µm long) microorganisms (potentially heterotrophic bacteria) attached to the cell walls of algae (**Figure 7C**). Beneath the algae zone, we found aggregates of cyanobacteria (**Figure 7D**) and an accumulation of heterotrophic

bacteria (**Figure 7E**). EPSs with different micromorphologies surrounded all cells. The EPSs associated with algal cells has a honeycomb-like structure (**Figure 7B**) whereas the EPSs in closecontact to cyanobacteria had a slime-like structure (**Figure 7D**). Both types of EPSs were loosely bound to the cells and did not take the shape of these cells. Groups of heterotrophic bacteria were surrounded by capsular-like polysaccharides (**Figure 7E**) that are distinguished according to their thickness and consistency (Rossi and De Philippis, 2015).

The structure of the gypsum hypoendolithic habitat, which did not contain algae, was significantly different than that of the surface cryptoendolithic habitat (**Figure 8**). Sizeable pores (0.2–2 mm in diameter) with compacted walls were located above a compact gypsum bottom layer (**Figure 8A**). The structure of these pores, together with the columns formed by lenticular gypsum crystals and with the dissolution features at the edges of the anhedral gypsum crystals, similar to the "Grating Fabric" defined by Artieda (2013), revealed distinct gypsum dissolution and crystallization processes. The microbial colonization in form of cyanobacteria aggregates covered the walls of the large pores and, unlike in the cryptoendolithic habitat, did not adhere to sepiolite nodules. Many cyanobacteria cell appeared empty and were most likely dead cells (**Figure 8B**).

#### Micromorphological Characterization of Algae

Growth of algal cells was only observed with solid BG-11 mineral medium. The growth was very slow and the great

extracellular polysaccharides. (A–C) Fractured algae cells with ultrastructural elements such as: cell walls (cw), chloroplasts (ch), lipoidal filling (lg), and gypsum particles (Gy); detached cell wall, asterisk; rod shaped microorganisms (potentially heterotrophic bacteria) attached to algae, white arrows; arrow head in (B) point to EPSs. (D) Image showing cyanobacteria (Cy) aggregates and slime EPSs (arrow head). (E) Image showing bacteria (B) evolved by capsular EPSs (arrow head).

majority of algal cells died within a few weeks of their cultivation. The cultivated cells were globular, uninuclear, 28–50µm in diameter and with smooth cell walls (**Figure 9A**). The cells possessed a reticulated chloroplast with a complex structure and often with several off-center pyrenoids (**Figure 9B**). Mature cells accumulated large quantities of orange-colored pigments (**Figure 9C**). Asexual reproduction occurred via 12–16 autospores, liberated by gelatinization of the maternal cell wall (**Figure 9D**). The autospores had irregular shape and size of 14–20µm in diameter. They possessed a chloroplast divided into several polygonal pieces lacking pyrenoids (**Figure 9E**).

# Molecular Characterization of Endolithic Microbial Communities

Taxonomic assignment revealed that the dominant phyla in both samples were Cyanobacteria, representing between 67 and 83% of all the sequence reads, and Proteobacteria, and Actinobacteria (**Figure 10A**). The only other photosynthetic phylum, in addition to Cyanobacteria, was that of Chloroflexi with 1.2–1.6% of sequence reads (**Figure 10A**). All diversity metrics, at a maximum sequencing depth of 1662 sequences per sample, revealed a higher diversity in sample CL1, with 260 OTUs at the 0.03% cutoff, when compared with sample CL2, with 167 OTUs (see Supplementary Table S2). Across the two samples, cyanobacteria included members of the Chroococcales, Synechococcales, Pseudanabaenales, and Nostocales; the most numerous were Chroococcales (64%) (mostly of the Chroococcidiopsis genus), followed by Synechococcales (18%), and unclassified Cyanobacteria (14%) (**Figure 10B**). No

FIGURE 8 | (A) SEM-BSE image of the bottom part of the gypcrete colonized by hypoendolithic cyanobacteria (Cy) forming aggregates; pore, p; sepiolite nodule, Sp; arrows point to cyanobacteria aggregates; hypoendolithic cyanobacteria colonization zone, Cy; bottom layer, Btm; gypsum, Gy. (B) Details of cyanobacteria aggregate.

FIGURE 9 | Morphology of cultivated green alga. (A) Size variability of vegetative cells. (B) Mature cell with a complex chloroplast and with several pyrenoids marked by arrows (in A,B). (C) Accumulation of orange color pigments. (D) Mature autosporangium. (E) Irregularly shaped autospores shortly after their liberation. Scale bars = 20 µm.

chloroplast sequences from algae were detected indicating a considerable phylogenetic distance away from the Cyanobacteria or potential primer mismatches with the plastid's 16S rRNA genes. This was not the case in a previous study were we obtained amplicons from the 16S rRNA gene of the chloroplast from an algae related to the Dolichomastix genera (Mamiellaceae family) with this same primer set used in this study (Robinson et al., 2015). Despite multiple attempts at PCR amplification, using a range of algal primers, with and without nested primers, we were unable to obtain the 18S rRNA gene sequence of the alga that was visualized by microscopy and analyzed by Raman spectroscopy. Future studies, using metagenome data from the gypsum community, together with our cultivation effort, might provide the genetic characterization of this alga.

The most abundant OTUs from each of the two gypsum datasets were used to build a phylogenetic tree. Most sequences belong to Cyanobacteria with closest relatives from arid deserts around the world and from members of the Chroococcidiopsis genus (Supplementary Figure S3). Only two OTUs were from Alphaproteobacteria, also associated with reference sequences from extreme environments (Supplementary Figure S3).

# Raman Spectroscopy of Endolithic Microbial Communities

We obtained Raman spectra for both the cryptoendolithic and hypoendolithic microbial habitats within the gypsum sample (**Figures 11A–D**). Micro-Raman spot analysis of orange (**Figure 11A**) and green (**Figure 11B**) cells from the cryptoendolithic zone revealed spectral records dominated by carotenoid and chlorophyll signals. No phycobiliproteins were detected suggesting that both spectra belong to algal cells. The spectra showed slight differences in wavenumber position of the carotenoid bands between the green and orange cells, with the ν1(C = C) band position at 1521 and 1519 cm−1 (orange cells at 785 and 532 nm) and 1524 and 1518 cm−1 (green cells at 785 and 532 nm). These differences were reproducible for other

FIGURE 11 | Examples of Raman spectra of orange algal cell (A), green algal cell (B) and cyanobacterial cells from the hypoendolithic habitat containing scytonemin (C) and without scytonemin (D). Spectra were obtained by the 785 nm excitation source; c, carotenoid; cx, carotenoid-xanthophyll; chl, chlorophyll; p, phycobiliproteins; s, scytonemin; gyp, gypsum. Asterisks (\*) indicate important corroborative bands for discrimination between scytonemin and phycobiliproteins. (E), microscopy image of cryptoendolithic habitat colonized by orange and green algae; (F) Raman distribution map of carotenoid (white spots) of the marked square in (E). For examples of Raman spectra extracted directly from the image see Supplementary Figure S5.

areas of the substrate and were interpreted as a result of selective resonance Raman effect (see Supplementary Section S8.4 section). The Raman distribution map of carotenoids in algae from the cryptoendolithic habitat displayed a higher intensity signal just beneath the gypsum surface and a lower intensity signal close to green algae zone (**Figure 11F**). The method of Raman imaging proved to be extremely powerful for mapping the spatial distribution of strong Raman scatters in geobiological samples (Vítek et al., 2014). Raman signals for cyanobacteria were clearly recognized in the gypsum hypoendolithic habitat. In particular for scytonemin (**Figure 11C**), a UV-screening pigment characterized by the corroborative Raman bands at 1595, 1554 cm−1 , and 1173 cm−1 (Edwards et al., 2000). The feature around 1322 cm−1 is the combination of scytonemin and chlorophyll signals. While brown-green cell aggregates produced a scytonemin signal, no scytonemin was detected within the neighboring blue-green cells (**Figure 11D**). Here, Raman features of carotenoids (at 1518 and 1156 cm−1 and Supplementary Figure S6 for Raman identification of carotenoids) and chlorophyll (1325, 915, and 743 cm−1 ) were detected. Additionally, strong bands for phycobiliproteins also revealed the presence of cyanobacteria. These accessory light harvesting pigments were characterized by the bands at 1627, 1582, 1370, 1283, 1237, 1108, 1050, 874, 820, 667, and 505 cm−1 (**Figure 11D**).

are observed. Carotenoids and scytonemin absorption ranges are indicated by double arrowed lines.

# Characterization of Metabolite Pigments by Spectrophotometry

The absorbance spectrum of pigment extracts obtained from different endolithic colonization zones are shown in **Figure 12**. Extracts from the orange algal cryptoendolithic zone showed high absorption values for the carotenoids absorption range. In contrast, extracts from the green cryptoendolithic colonization zone had lower absorption peaks for carotenoids and a higher chlorophyll absorption peak at 663 nm. This may indicate that the green zone also contains cyanobacteria together with green algae of lower carotenoid content. Extracts from the hypoendolithic zone showed strong absorption in the range of UVA corresponding to scytonemin absorption peaks at 315– 400 nm (Proteau et al., 1993). Semiquantitative determination, based on the analysis of absorbance spectrum from extracted pigments, showed that the gypsum orange zone had twice the content of carotenoids and less of both chlorophylls a and b than the green zone (Supplementary Table S2).

# Discussion

# How Extreme is the Environment at the Colonized Gypsum Site?

Most of the Atacama Desert is characterized by persistent extremely dry conditions and scarce precipitations in its hyperarid area (Houston and Hartley, 2003; Hartley et al., 2005; Houston, 2006a). Based on an evapotranspiration rate of 2920 mm y−1 for the nearby Salar de Atacama basin (Houston, 2006a), and yearly precipitations of 27.1 mm (this study), the aridity index for the Cordon de Lila sampling site is as low as 0.0093, presenting exceptionally hyperarid conditions. The extreme aridity of this location is also expressed by long timespans between scarce precipitations and also extremely low relative humidity. The mean annual RH value of 16.5% at Cordon de Lila is even lower than the value of 17.3% reported by Azua-Bustos et al. (2015) and recently reported as the driest site of the hyperarid Atacama Desert. Chile's northern Atacama Desert is also known as one of the places on Earth with the highest surface UV radiation (Cordero et al., 2014). The most extreme values for solar total and UV irradiances (1528 Wm−2 ) were reported for "Puna de Atacama" by Piacentini et al. (2003) and by Cabrol et al. (2014) (43.3 UV index), in a region of high altitude, with mostly cloudless conditions, relatively low total ozone column, and extremely low RH values. Recently, Rondanelli et al. (2015) show that the location of Earth's surface solar radiation maximum is in the Atacama Desert in the region between 24◦ and 25◦ S along 69◦W, over what is known as the Domeyko Cordillera. While solar radiation is essential for phototrophic life, extremely high levels of UV radiation is harmful to cellular components (Jeffrey et al., 1996; Phoenix et al., 2006), and excessive PAR results in photo-inhibition by damaging the reaction center of photosystem II in photosynthetic microorganisms (Häder, 1986). By combining both hyperaridity and extreme solar irradiance, the Cordon de Lila site with endolithically colonized gypsum studied here holds world records for environmental extremes, leading to the virtual absence of any life form on its rocks and soil surfaces.

# How Extreme is the Environment within the Gypcrete Endolithic Microhabitat?

The effectiveness of gypsum in attenuating solar irradiance has been evaluated by several studies. Amaral et al. (2007) reported UV transmission of less than 20% through a 0.5 mmthick gypsum layer. Effective attenuation of light, leading to transmission of 0.23% to 0.5% (from 300 to 680 nm) of incident light through a 0.9 mm thick Arctic gypsum crust, was reported by Cockell et al. (2010). The transmission of only 0.005% of the UVB and 0.05% of UVA, and 1% of PAR, were reported in a 1.2 mm thick Antarctic endolithically colonized gypsum by Hughes and Lawley (2003). In contrast, a 1 mm thickness of Arctic selenite exhibited only 25% reduction in ambient UV exposure (Parnell et al., 2004), underlying the importance of the micromorphology of gypsum crystals and the compactness of the gypsum layer in solar irradiance attenuation. Our data for the gypsum cryptoendolithic colonization zone, showing no UV transmission and 0.1–1% of PAR transmission, are consistent with previously reported studies. However, when taking into account the extreme solar irradiance at the sampling site and the small zenith angle, between 5◦ and 25◦ , experienced at the site when rain water is potentially available (July-March), it is likely that the cumulative dose of PAR within the endolithic habitat is extremely high.

Water is an essential element for the successful colonization of the gypsum deposits. The presence of cracks, fractures and pores in the gypsum substrate allow rainfall water and nutrients to penetrate into the gypcrete. Using EC measurements, we estimated that liquid water was present at the rock surface for less than 10 days per year. However, within the gypsum interior, the moisture regime could be significantly higher. We have observed compact structures and low porosity at the surface and bottom layers of the gypsum deposits that might significantly decrease water evaporation rates. We also identified gypsum crystal dissolution and crystallization patterns mostly at the bottom layer, indicating the presence of liquid water, within the rock, for some periods of time. These observations lead to what we call the "eggshell" shield effect in which the gypsum interior is protected against rapid water evaporation rates. Moreover, the water can be efficiently trapped and retained inside the gypsum deposits by the sepiolite nodules randomly distributed in the substrate. Sepiolite is a fibrous sheet of magnesium silicate clay with unique physicochemical and rheological properties due to the presence of zeolitic channels in its structure (Leguey et al., 2014). Its nano-porosity induces efficient water absorption and retention up to 0.3 g g−1 dry weight (Caturla et al., 1999). We suggest that the sepiolite water retention capacity might increase significantly the fraction of water available to endolithic microorganisms.

# Sub-micro Endolithic Habitats and their Colonizers

Our molecular characterization of the gypsum endolithic microbial community revealed a community dominated by cyanobacteria, mostly of the Chroococcidiopsis genus. The difference in community structure and composition we observed between samples suggest a high level of heterogeneity in the substrate's architecture, providing sub-micro habitats for colonizing microorganisms. A model of the gypsum deposit "architecture" with its cryptoendolithic and hypoendolithic habitats, gypsum layers, and sepiolite nodules, together with the role played by these structures and by the prevalence of gaseous or liquid water is illustrated in **Figure 13**.

In this model, we hypothesize that the water prevalence states and therefore the moisture regime inside the gypcrete is dependent on the depth, fabric, and others accompany to gypsum minerals. Following rainfall, water fills the fissures of the rock, forming percolation pipes and conducting the water toward the gypsum's bottom layer. There, the water, or rather a calcium sulfate solution, accumulates for longer times, filling the large pores or "micro-caves" we observed. The processes of gypsum dissolution-recrystallization and water evaporation convert the bottom layer of the gypcrete into a compact, thick layer with low permeability and potentially higher UV transmission rates (Parnell et al., 2004). Liquid water is also retained within sepiolite nodules and pore spaces inside the gypcrete, providing sub-micro environments with high moisture until the gypsum is completely dried out.

The first adaptation strategy for microorganisms in hyperarid deserts is the avoidance of extreme external conditions by escaping into the endolithic habitat (Karsten and Holzinger, 2014). The second step is to follow the water, which is essential for the metabolic activity of the primary producers in the community, followed by the level of PAR. This explains the well-defined colonization patterns we observed for the algae and cyanobacteria. We propose that moisture is the key abiotic factor driving the spatial distribution of chlorophototrophs in the gypsum sub-micro habitats, providing there are sufficient PAR levels. The cryptoendolithic habitat is an area with scarce liquid

deposits from the Atacama Desert.

water with the exception of sepiolite moisture islands and pores, with high RH values, among gypsum crystals. Our results show that this habitat is dominated by a layer of green/orange algae and a subjacent layer of cyanobacteria attached to sepiolite nodules. Cultivation of the endolithic green algae showed that the orangeand green-colored layers were formed by the same organism. This alga very probably belongs to the class Chlorophyceae (Chlorophyta), known to contain several desert soil algae (Lewis and Lewis, 2005). Morphologically, the gypsum alga resembles the genera Spongiochloris, Neochloris, and Deasonia. However, algae from natural gypsum cryptoendolithic habitat exhibit a distinct suite of morphological features that is unique among all described algal genera. Therefore, we conclude that the alga observed in this study belongs to a new, yet undescribed genus of Chlorophycean green algae.

Previous works indicated that, in arid environments, high relative humidity was enough to induce metabolic activity in algae. For example, photobionts from cryptoendolithic lichens inhabiting sandstones of the Dry Valleys of Antarctica begin to photosynthesize at RH = 70% (Palmer and Friedmann, 1990a) and cultured algae (Trebouxia sp.) are capable of photosynthesis at a RH of 80% (Palmer and Friedmann, 1990b). In contrast, endolithic cyanobacteria inhabiting sandstones in the Negev Desert require RH > 90% for photosynthesis (Palmer and Friedmann, 1990a). These findings may explain the close contact of cyanobacteria with sepiolite nodules within the gypsum. However, the hypoendolithic habitat in the gypsum deposits is an area with periodical occurrence of liquid water after rainfalls hence it is colonized by cyanobacteria colonies covering the walls of micro-caves. The hypoendolithic habitat was first described within the gypsum crusts collected in a less extreme (from the point of view of high night RH values) zone of the Atacama Desert (Wierzchos et al., 2011; Vítek et al., 2013). However, these authors reported the presence of algal Trebouxialike cells aggregates closely associated with fungal hyphae as the main hypoendolithic colonizers of the bottom gypsum crusts, suggesting scarce liquid water within that hypoendolithic habitat.

Light penetration has been considered a primary control for the location and extent of colonization by endolithic phototrophic microorganisms (Walker and Pace, 2007). Endolithic phototrophs from Antarctica (gypsum and sandstone), and other regions around the world, were found to receive as little as 0.01% of incident light (Nienow et al., 1988; Matthes et al., 2001; Hughes and Lawley, 2003). In our study, we did not observed chlorophototrophs below 0.1% of incident PAR in the cryptoendolithic habitat. We suggest that the additional hyperaridity found at the Cordon de Lila sampling site could potential hinder the colonization of deeper substrate layers where PAR levels are below 0.1%. However, despite relatively low values of transmitted PAR, absolutely PAR values were as high as 30µmol photons m−2 s −1 , potentially inducing chlorophototrophs to produce abundant secondary metabolites such as carotenoids and scytonemin. High rates of solar irradiance increase considerably the probability of photoinhibition and photooxidative lethal damage to phototrophs (Solovchenko and Merzlyak, 2008). Photoinhibition has been defined as the inhibition of photosynthesis caused by excessive radiance; it may damage the photosynthetic apparatus, causing the (photo-) destruction of the photosynthesizing pigments (Long et al., 1994; Alves et al., 2002). It is therefore likely that the algae and cyanobacteria colonizing the gypsum endolithic habitat have evolved self-protective strategies such as the efficient production of secondary metabolites shielding them against UV and PAR excessive irradiance.

# Solar Irradiation Protective Pigments Nature and their Spatial Distribution

Carotenoids are the accessory pigments of the photosynthetic apparatus, which participate in both light-harvesting and photoprotective functions (Nguyen et al., 2013). Some shortchained carotenoids can potentially act directly as UV-screening pigments but are not favored by organisms, whereas ubiquitous carotenoids with nine or more conjugated double-bonds absorb more in the light visible region (Cockell and Knowland, 1999). These long-chain carotenoids are able to cope with toxic reactive oxygen species that can form during photosynthesis, or as the result of UV-irradiation, and which would destroy the biomolecules forming the photosynthetic apparatus (Asada, 1994; Telfer et al., 2008).

From our Raman spectroscopic results, we observe a clear enhancement of carotenoid/chlorophyll intensity ratio of Raman signal in the case of the orange-color algal cells compared to the green-color algal cells. All the band positions may be attributed to a standard C-40 carotenoids, with bands above 1520 cm−1 probably due to xanthophyll compounds such as zeaxanthin and lutein, and the band at 1514 cm−1 assignable to β-carotene. However, based just on the Raman spectra alone the exact and unequivocal identification of a particular carotenoid in these cells was not possible, as the ν1 band position can be identical for several carotenoids or different for the same carotenoid in different biological matter (de Oliveira et al., 2010). We interpret the differences in relative carotenoid content (increased within orange algal zone) and carotenoid composition as a response of algae to light-driven stress enhanced with the proximity to the surface. The relative intensity of chlorophyll Raman spectral features differs significantly between green and orange algal cells. The clear characteristic signals of chlorophyll were detected within the green cells using 785 nm laser for excitation, represented by the bands at 1326, 918, and 745 cm−1 . On the other hand, spectra from the orange cells revealed none or very weak chlorophyll features. The very weak signal of chlorophyll in the orange cells was the result of decreased chloroplasts content in these cells, which was clearly confirmed by microscopy observations.

The capability to synthesize very high amounts of a complex mixture of secondary carotenoids under environmental stresses is widely spread over the green microalgae, as a defense mechanism against environmental injuries. Secondary carotenoid formation enhances the cell resistance to oxidative stress generated by unfavorable environmental conditions including excess light, UV and PAR irradiation, and nutrition stress and, therefore, confers a higher survival capacity to the cells (Hanagata and Dubinsky, 1999; Orosa et al., 2000, 2001). High light levels and salinity were reported as main stress drivers enhancing carotenogenesis in the aerial microalga Scenedemus sp. (Chlorophyceae) leading to the production of xanthophylls as fatty acid esters in the algal cells (Aburai et al., 2015). Moreover carotenoids have been identified as photoprotection agents against photoinhibition and molecular damage of the photosynthetic apparatus (Oren et al., 1995).

Our analytical results obtained by spectrophotometry and Raman spectroscopy demonstrate enhancement in carotenoids and depletion in chlorophyll photosynthetic pigment within orange cells lying close to the gypsum surface. This progressive carotenogenesis process in the case of the cryptoendolithic algae from the Atacama's gypsum deposits is driven by higher PAR light levels toward the gypsum surface, the oligotrophic environment, and the very large periods (up to 9 months) of dehydration followed by desiccation. It was also described that carotenogenesis processes might lead certain algal cells to turn from a green to an orange/red-color (e.g., "blood alga" Haematococcus pluvialis, snow alga Chlamydomonas nivalis or Scenedesmus komareki) due to the accumulation of large amount of carotenoids within the cell (Hanagata and Dubinsky, 1999 see also references therein).

Carotenogenesis generally accompanies larger metabolic changes and morphological modifications, i.e., progressive accumulation of lipoidal globules and/or disappearance of certain cell organelles such as chloroplasts (Hanagata and Dubinsky, 1999). A similar process has been observed in the case of the cryptoendolithic algae, where lipoidal globules that were first formed at the periphery of the cell, progressively propagated toward its inside, eventually filling most of it. The chloroplast was single and parietal in the green cells. As the cells turned orange, the chloroplast divided into several small lobes and was pushed toward the interior of the cell by the accumulating lipoidal globules. However, our knowledge of the cellular localization of secondary carotenoids and the change of cellular ultrastructure associated with their accumulation in other carotenoid-accumulating algae is limited and will be the topic of future studies.

The presence of scytonemin, a well-known UV-screening pigment, in the gypsum hypoendolithic habitat was confirmed by Raman spectroscopy and spectrophotometric analyses. Scytonemin is a typical cyanobacterial metabolite that is extracellularly located in the cyanobacteria exopolysaccharide sheath around cells and cell aggregates, as recently shown for endolithic cyanobacteria colonizing Atacama's halite (Vítek et al., 2014). Scytonemin is known to have photoprotective properties against harmful UV radiation (Garcia-Pichel and Castenholz, 1991; Cockell and Knowland, 1999; Vítek et al., 2014). Periodic desiccation of cyanobacteria in combination with high UV fluxes may further induce scytonemin biosynthesis (Fleming and Castenholz, 2007). This raises the question of how can UV irradiate hypoendolithic cyanobacteria? First, it should be mentioned that hypoendolithic colonization has only been found in these specific gypsum deposits, which were detached from the coarse grain soil underneath, suggesting that diffused solar light could illuminate the gypcrete underneath. In addition, as mentioned above, the bottom gypsum layer shows a very compact structure potentially permitting some transmission of UV (Parnell et al., 2004) from intense solar irradiance.

# EPSs and their Role in Reducing Water Stress and UV Radiation

Our results show abundant presence of exopolysaccharides with different micromorphology surrounding algae, cyanobacteria, and heterotrophic bacteria single cells and colonies. EPSs are an important class of biopolymers, which play key protective and structural roles. Under extremely low water availability, large periods of dehydration, and excessive solar irradiance, the protection against water stress and UV and PAR radiation is one of the main roles of EPSs in constrained environments (Rossi and De Philippis, 2015). Although the role of EPSs in reducing water stress has not been fully clarified, they are reportedly involved in maintaining hydration because of their hydrophilic/hydrophobic characteristics (Makhalanyane et al., 2015). Hygroscopic EPSs act as a water reservoir after wetting events (Gorbushina, 2007) and as a possible mechanism for absorption and retention of water vapor, as reported for hypolithic microbial communities from the Dry Valleys of Antarctica (de los Ríos et al., 2014b). In addition, it was described that EPSs may contain UV-absorbing compounds, such as scytonemin (Vítek et al., 2014).

# Architecture of a Lithic Substrate

As observed with microscopy techniques, internal structural elements such as porosity, pore-size distribution, presence of large pores and cavities, light transparency, and light scattering properties, dissolution and crystallization features, and sepiolite nodules distribution varied significantly within various locations of the gypsum substrate. All together these structural, physical, chemical and mineral properties give rise to a new understanding of the features and functions relevant to the rock bioreceptive characteristics. We suggest using the term of "rock architecture" instead of "rock structure" to emphasize the functional role of the rock interior. As such, this new concept of the architecture of a lithic substrate encompasses all the internal structures of a rock that are essential as a habitat for microbial life. It is about perceiving the rock interior from the existence of porous spaces of different sizes and shapes, interconnected or not; the solid structures that divide and support these spaces, and the minerals and salts that can be transformed. All these components and elements are interrelated and influence one another, thus fulfilling a requisite: they shape a suitable architecture to hold microbial life. The architecture of habitable rocks provides resources (water, light, and nutrients, above all) and guarantees effective protection from excessive evapotranspiration, thus assuring efficient gas exchange. It also provides a long-term stable environment for life. Considering the architecture of a rock provides an integrated view of its potential habitability for endolithic microbial communities. All porous rocks have a structure, yet very few show such a suitable architecture for endolithic microbial colonization, even under extreme environmental conditions, as the Atacama's gypsum do.

# Concluding Remarks

Here we report that architectural features of translucent gypsum rocks with sepiolite inclusions provide increased water availability and attenuation of harmful UV and PAR radiation to endolithic microbial communities. Intense PAR radiation induced chlorophototrophic microorganisms to form a unique layering of algae and cyanobacteria within the cryptoendolithic habitat and the synthesis and accumulation of carotenoids in the upper algae zone. Scytonemin pigment is produced by the hypoendolithic cyanobacteria to avoid photoinhibition and photooxidative damages. We found that the spatial distribution of photosynthetic pigments within the gypsum rock such as carotenoids, chlorophylls and phycobiliproteins is linked to different colonization zones and types of microorganisms. The highest concentration of carotenoids was detected in orange-colored algal cells just beneath the gypsum surface, suggesting that the green microalgae use this shielding as a survival strategy against the intensive solar irradiance. This in turn allows green microalgae to colonize a slightly deeper layer in the rock until maturation and transformation into orange-colored algae with high carotenoid content. The same carotenoid-like "umbrella" is also protecting cyanobacteria cells underneath the layer of green microalgae cells. The detection of scytonemin in the hypoendolithic habitat colonized by cyanobacteria can be interpreted as another adaptation strategy

# References

against excess of UV and PAR in this unique conformation. The highly specialized microbial assemblages we describe here can be use as model systems to further our understanding of photoprotective mechanisms under extremely unfavorable conditions.

# Acknowledgments

CA, JD, OA, AD, and JW are thankful for financial support by CGL2013-42509P grant from MINECO, Spain. JD acknowledges funding from the NASA Exobiology Program (Grant EXOB08– 0033) and from the National Science Foundation (Grant NSF-0918907), PV acknowledges funding from the Czech Science Foundation, project no. P210/12/P330 and ENVIMET project CZ.1.07/2.3.00/20.0246 and AD acknowledges funding from the NASA Exobiology Program (Grant NNX12AD61G) and the NASA Astrobiology Institute (NAI Grant NNX15BB01A to the SETI Institute). The authors would like to thank to the MNCN - CSIC Microscopy Service staff and A. Gonzalez (CIC – Univ. Granada) for technical assistance.

# Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fmicb. 2015.00934

halite evaporites of the Atacama Desert. Inter. Microbiol. 13, 79–89. doi: 10.2436/20.1501.01.113

**Conflict of Interest Statement:** 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.

Copyright © 2015 Wierzchos, DiRuggiero, Vítek, Artieda, Souza-Egipsy, Škaloud, Tisza, Davila, Vílchez, Garbayo and Ascaso. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Namib Desert dune/interdune transects exhibit habitat-specific edaphic bacterial communities

Sandra Ronca1 , Jean-Baptiste Ramond1 , Brian E. Jones 1, 2, Mary Seely 3, 4 and Don A. Cowan1 \*

1 Department of Genetics, Centre for Microbial Ecology and Genomics, University of Pretoria, Pretoria, South Africa, 2 DuPont Industrial Biosciences, Leiden, Netherlands, 3 Gobabeb Research and Training Centre, Walvis Bay, Namibia, 4 School of Animal, Plant and Environmental Sciences, University of the Witwatersrand, Johannesburg, South Africa

#### *Edited by:*

Paul Bodelier, Netherlands Institute of Ecology (NIOO-KNAW), Netherlands

#### *Reviewed by:*

Nathan Basiliko, Laurentian University, Canada Roey Angel, University of Vienna, Austria Sarah Elizabeth Evans, Kellogg Biological Station Michigan State University, USA

#### *\*Correspondence:*

Don A. Cowan, Centre for Microbial Ecology and Genomics, University of Pretoria, Natural Sciences Building II, Hatfield 00028, South Africa<EMAIL_ADDRESS>
#### *Specialty section:*

This article was submitted to Terrestrial Microbiology, a section of the journal Frontiers in Microbiology

*Received:* 29 April 2015 *Accepted:* 03 August 2015 *Published:* 04 September 2015

#### *Citation:*

Ronca S, Ramond J-B, Jones BE, Seely M and Cowan DA (2015) Namib Desert dune/interdune transects exhibit habitat-specific edaphic bacterial communities. Front. Microbiol. 6:845. doi: 10.3389/fmicb.2015.00845 The sand dunes and inter-dune zones of the hyper-arid central Namib Desert represent heterogeneous soil habitats. As little is known about their indigenous edaphic bacterial communities, we aimed to evaluate their diversity and factors of assembly and hypothesized that soil physicochemistry gradients would strongly shape dune/interdune communities. We sampled a total of 125 samples from 5 parallel dune/interdune transects and characterized 21 physico-chemical edaphic parameters coupled with 16S rRNA gene bacterial community fingerprinting using T-RFLP and 454 pyrosequencing. Multivariate analyses of T-RFLP data showed significantly different bacterial communities, related to physico-chemical gradients, in four distinct dune habitats: the dune top, slope, base and interdune zones. Pyrosequencing of 16S rRNA gene amplicon sets showed that each dune zone presented a unique phylogenetic profile, suggesting a high degree of environmental selection. The combined results strongly infer that habitat filtering is an important factor shaping Namib Desert dune bacterial communities, with habitat stability, soil texture and mineral and nutrient contents being the main environmental drivers of bacterial community structures.

Keywords: deterministic assembly, open soil, habitat filtration, dune transects, environmental gradient, arid ecosystems

# Introduction

Despite extreme surface conditions, notably (hyper-) aridity, wide daily temperature fluctuations, high UV radiation and oligotrophy, edaphic microbial communities have been shown to flourish in desert soils (Makhalanyane et al., 2015). The Namib Desert of south-western Africa is considered to be one of the oldest desert regions on Earth (Eckardt and Spiro, 1999). It is characterized by a wide range of different soil environments including gravel plains, sand dunes, inselbergs, escarpments, river beds, salt pans and playas (Seely, 2012; Eckardt et al., 2013a). In the central Namib Desert, studies have notably determined that local soil physicochemical conditions and climate played a significant role in the assembly of edaphic and hypolithic communities were (e.g., Stomeo et al., 2013; Warren-Rhodes et al., 2013; Ramond et al., 2014; Gombeer et al., 2015).

However, and despite representing over 41% of the total desert land surface (Seely, 2012), the detailed microbial ecology (as determined by modern molecular tools) of the Namib Sand Sea has to date poorly been assessed and essentially focused on fungal communities (Jacobson, 1997; Jacobson and Jacobson, 1998; Jacobson et al., 2015). The dispersal and colonization of mycorrhizal fungal communities were shown to depend on sand stability and moisture availability, respectively (Jacobson, 1997). And, moisture-activated decomposing fungi were observed on the litter of the perennial Namib dune grass Stipagrostis sabulicola and in dune sands (Jacobson and Jacobson, 1998; Jacobson et al., 2015). Recently, a T-RFLP fingerprinting-based study compared the edaphic bacterial community assemblages from multiple soil surface geologies and lithologies of the central Namib Desert. It showed that (i) interdune bacterial communities were significantly different than those of the gravel plains and riverbeds and that (ii) the soil physicochemistry and lithology of the different interdunes sampled were important in structuring their indigenous communities (Gombeer et al., 2015).

Dune morphology and dynamics have in contrast extensively been studied in the Namib Desert. Three geographically distinct dune morphology patterns have been characterized: transverse dunes to the west (i.e., toward the coast), linear and complex linear dunes in the center and star dunes in the east (i.e., inland) (Livingstone, 2013). Similarly, an east/west sand color gradient, from yellowish brown (west) to much redder (east), has been described and attributed to variations in iron oxide and clay mineral concentrations (Walden and White, 1997; Livingstone, 2013). Within dune slopes, a complex moisture gradient has been shown to influence the number of species perennial grasses as well as their growth strategies and distribution (Yeaton, 1988; Seely, 1990). Dunes have also been described to be highly dynamic environments, as they are constantly modified by surface shear forces and mechanical stresses. For example, it has been reported annual transverse dune crest migrations varying from 4 to 56 m (Livingstone, 2013). Wind has also been shown to alter the top of wind-exposed slopes (or stoss slopes) and sand accumulation to dominate along lee slopes (or downwind slopes) (Lancaster, 1985; Eckardt et al., 2013a). The interactions between wind erosion and sediment deposition lead to discernible patterns in grain size and sorting gradients both on individual dune slopes and between dune and interdune areas.

In such a dynamic and heterogeneous habitat, dune microbial communities have thus to cope with numerous stresses, including environmental physical instability, fluctuating soil physicochemical properties, low water availability, high temperatures, oligotrophy and alkaline pH (typically around 8) (Makhalanyane et al., 2015); all of which have been shown to affect both community diversity and function (Andrew et al., 2012; Kuske et al., 2012; Yu and Steinberger, 2012a; Makhalanyane et al., 2015).

Using the Kahani dune system in the Namib Sand Sea as a model dune environment, we aimed to determine drivers shaping edaphic bacterial community diversity in four dune biotopes (dune top, slope, base, and interdune). We hypothesized that in such a dynamic environment, strong soil physicochemical gradients would structure bacterial communities; i.e., that deterministic processes would lead to dune habitat-specific communities. Our related null hypothesis stated that high dispersal rates, due to the important dynamicity of dune ecosystems, would lead to homogenized communities (Leibold and Norberg, 2004). To test these hypotheses, we collected 125 soil samples from 5 parallel dune/interdune transects in a single Namib Desert dune system; i.e., from the summit of a western dune to the crest of an eastern dune. For each sample, we recorded 21 edaphic physico-chemical parameters and characterized the bacterial community structure by 16S rRNA gene T-RFLP analyses. We also performed 16S rRNA gene pyrosequencing on representative samples of 7 putative dune habitats found across the dune/interdune transect, namely the east and west dune tops, slopes bases and the interdune, to investigate their taxonomic composition and assess if they presented different assemblages.

# Materials and Methods

# Study Site and Soil Sample Collection

Soil samples were collected in April 2013 in the Kahani dune system (S 23◦ 35.830′ /E 15◦ 01.800′ ), a typical complex linear dune in the Namib Desert Sand Sea (**Figure 1**). 25 soil samples were collected from each of five parallel (40 m apart) 2 km long dune/interdune linear transects from the west (stoss slope) to the opposite east dune (lee slope) resulting in a total of 125 soil samples (**Figure 1**). Each transect was divided into 4 morphologically distinct zones: dune Top, Slope, Base and Interdune. Sampling points were regularly distributed; every 40 m from the west dune Top to its base, every 200 m in the Interdune area and every 20 m on the eastern dune up to the dune Top, as the west facing dune was typically shorter and steeper than the east facing dune (**Figure 1**). Surface soil samples (0- 5 cm depth) were collected (avoiding grasses on the dune Slope to minimize any rhizosphere bias), and each individual sample comprised the pooling of four pseudo-replicate soils in a 1 m2 area. Samples were stored at 4◦C for soil chemistry analysis, at −20◦C for molecular analysis and processed immediately for fluorescein di-acetate (FDA) assays.

# Soil Physico-chemistry Analyses

All analyses were performed at the Soil Science Laboratory of the University of Pretoria. Prior to analyses of 21 edaphic parameters, the soil samples were sieved (2 mm mesh) and dried overnight at 37◦C. pH was measured in 10 g soil slurries (1:2.5 soil/deionized water) with a pH meter (Crison basic 20, Barcelona, Spain) (Eckert and Sims, 1995). Soil nitrate (NO− 3 -N) and ammonium (NH+ 4 -N) were determined by extraction (2M KCl) and steam distillation with subsequent titration performed as described by Keeney and Nelson (1982). Soil total carbon percentage (C) was measured by oxidizing the organic material with chromic acid and titrating the excess dichromate (Walkley, 1935; Nelson, 1982). Cation exchange capacity (CEC) was determined using an ammonium acetate solution as extractant for the exchangeable and water-soluble cations with steam distillation to separate the ammonia (Rhoades, 1982). Ammonium acetate extraction was also used to measure salt concentrations (Na, K, Mg, and Ca) using inductively coupled plasma atomic emission spectroscopy (ICP-OES) (Spectro Genesis, Spectro Analytical Instruments GmbH, Germany). The Bray-1 method (Bray and Kurtz, 1945) was used to quantify extractable phosphorous (P). Metals (Fe, Al, and Mn) were extracted from 10 g of soil using EDTA buffer and quantified by ICP-OES after filtration through a of the barren Interdune.

0.45µm Millipore filter (EMD Millipore Corporation, Billerica, MA, USA).

Particle size distributions were determined using the ASTM protocol D422-63 on 50 g of soil sample (ASTM, 1985). The hydrometer method (Bouyoucos, 1962) was used to measure the size distribution of the silt and clay fractions obtained from the ASTM protocol. Granulometry analysis contributed to the formation of 8 variables: very coarse sand (VCS), coarse sand (CS), medium sand (MS), fine sand (FS), very fine sand (VFS) and the portion of the silt plus clay (SC), silt (ST), and clay (CY).

# Fluorescein Di-acetate Degradation Assay

FDA hydrolysis, as an indicator of general bacterial metabolic activity, was measured using the protocol of Green et al. (2006). 0.5 g of soil was mixed with 12.5 mL of sterile phosphate buffered saline (PBS, pH 7.4) and 0.25 mL of 4.9 mM FDA dissolved in acetone and subsequently incubated at 28◦C for 1 h under constant agitation. FDA hydrolysis was stopped by mixing 40µL of acetone with 1 mL of soil slurry. After centrifugation (8800 × g, 5 min.) fluorescence was measured with a portable fluorometer (Quantifluor, Promega, Madison, USA), calibrated against a standard curve. FDA activity was expressed as µg fluorescein per g soil per hour (µg fluorescein g−1 soil h−1 ).

# Metagenomic DNA Extraction, PCR Amplification and T-RFLP

Metagenomic DNA was extracted from 0.5 g of soil using the PowerSoil™ DNA Isolation Kit (MoBio, West Carlsbad, CA, USA) following the manufacturer's instructions.

PCR reactions were carried out in a Thermo Cycler (BioRad) in 50µL final volume reactions containing 1X PCR Buffer A, 1X Enhancer, 0.2 mM each dNTP, 0.5µM of each primer (E9F: 5′ -GAGTTTGATCCTGGCTCAG-3′ /U1510R: 5′ - GGTTACCTTGTTACGACTT-3′ ) (Reysenbach and Pace, 1995; Hansen et al., 1998), 4% DMSO, 0.5 U of KAPA2G Robust DNA polymerase (KAPA, USA) and 10 ng of template DNA. Thermal cycling conditions were as follows: 3 min denaturation at 95◦C followed by 30 cycles with denaturation at 95◦C for 30 s, annealing at 55◦C for 30 s, and elongation at 72◦C for 90 s with a final elongation at 72◦C for 10 min. This reaction yielded 1.5 kb amplicons.

For T-RFLP analysis, the forward E9F primer was FAMlabeled and the PCR products were purified using the NucleoSpin Kit™ (Macherey-Nagel, Germany) and digested overnight using HaeIII (Fermentas, USA). After a second purification, electrophoretic separation of the terminal-restriction fragments (T-RFs) was performed using an ABI3130XL sequencer (Applied Biosystems, USA). The retrieved T-RFLP profiles were analyzed using Peak Scanner 1.0 (Applied Biosystems; http://www. appliedbiosystems.com). True peaks and fragments of similar size were identified and binned using the software R and Perl (Abdo et al., 2006). The term OTU (Operational Taxonomic Unit) is used to refer to individual terminal restriction fragments (T-RF) in T-RFLP patterns, with recognition that each OTU may comprise more than one distinct bacterial ribotype.

# Statistical Analyses

Principal Component Analyses (PCA) was performed (prcomp R function) on standardized soil chemistry variables and granulometric parameters (Oksanen, 2013). T-RFLP data were Hellinger-transformed (Legendre and Gallagher, 2001) and used to calculate Bray-Curtis dissimilarity matrices, which were further visualized using 3D non-metric multidimensional scaling (3D-nMDS) (Primer 6 software; Primer-E Ltd, UK). Permutational multivariate analysis of variance (PERMANOVA), function Adonis (vegan package for R), was used to test for significant differences between sample groups along the dune/interdune transects. Distance based redundancy analysis (dbRDA) was used to determine correlations between bacterial community structures and habitat parameters (e.g., soil chemistry and granulometry) (vegan package for R). Negative eigenvalues were transformed by taking the square root of dissimilarities.

An analysis of multivariate homogeneity of group dispersions (variances) (betadisper function in vegan package for R) was used to test if one or more of the dune zones was significantly more variable than the others. The test was performed on standardized soil physico-chemical properties and bacterial community structures (relative abundance of OTUs).

Edaphic variables significantly correlated, or "redundant" (based on PCA analysis), were grouped for dbRDA analysis in order to decrease co-linearity and noise (**Table 1**) (Legendre and Legendre, 2012), yielding 15 independent variables: pH, NO− 3 - N, NH+ 4 -N, C, CEC, Na, K, P, Mg (also representing Ca), Fe (also representing Al and Mn), very coarse sand (VCS), coarse

All values are given as the mean of ± standard error. In bold are indicated the variables grouped for the RDA analysis presented in *Figure 7*.

sand (CS), medium sand (MS), fine sand (FS) and very fine sand (VFS) (representing silt plus clay [SC], silt [ST] and clay [CY]). Tukey's HSD (honest significant difference) test was used to test for differences in averaged bacterial activity between dune areas.

#### 454 Pyrosequencing

Extracted metagenomic DNA samples from 7 dune zones (Tops E+W, Slopes E+W, Bases E+W and Interdune) were pooled. Bacterial 16S rRNA gene amplicons were generated with the primer set 27Fmod (5′ -AGRGTTTGATCMTGGCTCAG-3 ′ )/519Rmodbio (5′ -GTNTTACNGCGGCKGCTG-3′ ) and sequenced using a Roche 454 FLX titanium next-generation sequencer by Mr DNA Laboratories (Shallowater, Texas, USA). The sequences are available at the NCBI Sequence Read Archive under the accession number SRP059482.

#### Pyrosequencing Data Analyses

Analyses were carried out using the mothur software package (Schloss et al., 2009). Sequences with ambiguity, and/or with homopolymeric stretches longer than 8 bp, or shorter than 200 bp were removed with no barcode and primer mismatches accepted. This resulted in reads with sizes ranging from 232 to 314 bp. Chimeras were removed using Perseus (Quince et al., 2011). Sequences were aligned to the reference SILVA database and clustered into operational taxonomic units (OTUs) at the species (97% similarity) level using the average neighbor settings in mothur. Each OTU was assigned a taxonomic classification with reference to the Ribosomal Database Project (RDP; http:// rdp.cme.msu.edu/index.jsp). Collector's curves were produced using the Chao1 diversity index and the coverage using Good's coverage estimator. The community composition of the sand dune zones was compared with a UPGMA tree based on the Jaccard index. Indicator species analysis was conducted using the multipat function of the indicspecies package in R.

# Results

## Dune Soil Physico-chemistries

Principal Component Analysis of 21 physico-chemical parameters from all 125 independent soil samples (**Figure 2A**) along the five dune-interdune transects (**Figure 1**) shows Top and Slope samples clearly separated along the PC2 axis (**Figure 2A**), explaining 10.7% of the sample variation. Dune Base samples cluster between Slope and Interdune samples (**Figure 2A**). Top and Slope samples clearly separate from those of the Interdune along PC1, which explains the majority of physico-chemical soil variation (48.4%). The most significant result to emerge from our analysis is the observation that each of the four sites examined, namely [dune] Top, Slope, Base and Interdune, was distinct and significantly different (adonis R 2 = 0.477, p < 0.001) with respect to soil physico-chemistry. The test for multivariate homogeneity of group dispersions (betadisper) showed that variations in soil physico-chemical

properties within dune zones were not significantly different (p = 0.14). No statistically significant differences in physicochemical soil properties were found between dunes slopes exposed to different wind regimes (i.e., the lee or stoss slopes) (data not shown).

Soil texture analysis (**Figure 2B**) also shows separation of the four dune zones with an important fraction of fine sand in the Top zone and proportionally more coarse and medium sands in the Slope and Base areas (**Table 1**). The Interdune was characterized by a higher proportion of very fine sand, silt and clay than the other zones of the dune (**Table 1**). Dune Top samples had an elevated Na content and a higher pH (**Table 1**). Concentrations of minerals (K, Ca and Mg) and metals (Fe, Al) were lowest on dune Slopes (**Figure 2C**) and highest in the Interdune samples (**Table 1**) and inversely correlated with pH (**Figure 2C**). All other edaphic properties measured, including C, P, pH, CEC, NO− 3 -N, NH+ 4 -N, and VCS (**Table 1**), were only very weakly correlated with the PC1 axis (**Figure 2D**) and do not appear to contribute to the dune/interdune physico-chemical gradients. Data for the dune Base (**Table 1**) was generally intermediate between Slope and Interdune. To conclude, this comprehensive dune soil physico-chemistry analyses clearly defined four dune biotopes; i.e., Top, Slope, Base and Interdune.

# Dune Zone-specific Bacterial Community Structure and Diversity

16S rRNA gene T-RFLP analysis yielded between 12 and 105 OTUs per individual sample (α-diversity) for each of the 125 samples. The betadisper function for multivariate homogeneity indicated that variances in the relative abundance of OTUs within dune zones were not significant (p = 0.59). Dune Tops showed the lowest bacterial diversity (41 OTUs), while the Slope was the most diverse (80 OTUs) followed by the Base and Interdune with 62 and 53 OTUs, respectively. Only 16 OTUs (15%) were shared by the bacterial communities of the four dune biotopes (Figure S1A), with samples from the Slope showing the highest number of unique OTUs (16). The 3D-NMDS plot showed that dune bacterial communities cluster according to their soil biotope of origin; i.e., Top, Slope, Base and Interdune (**Figure 3**). This was confirmed by PERMANOVA analysis (adonis R2 = 0.31, p = 0.001) and consistent with the analysis of soil physico-chemistry (**Figure 2A**).

Pyrosequencing analysis of 16S rRNA gene PCR amplicons yielded a total of 5204 OTUs (defined at 97% sequence similarity, **Table 2**; Figure S2B), of which 2741 (47.3%) were singletons. Dune Tops (E+W) showed the highest number of unique OTUs (1154) followed by the dune Bases (1170), Slopes (971), and the Interdune (775) samples. Only 88 OTUs (2%) were shared

TABLE 2 | Distribution and diversity of OTUs (97% cut-off) in the seven dune zones studied.

E, East; W, West.

between the four dune biotopes, which is broadly consistent with the T-RFLP analysis (Figures S1A,B).

After resampling this dataset for sequence consistency (n = 4653 reads), collector's curves using the Chao1 index indicated that dune Tops (E+W) and Slopes (E+W) samples reached near saturation (Figure S2). UPGMA analysis shows that bacterial communities cluster strongly in relation to their dune biotope of origin (**Figure 4**), supporting the T-RFLP analysis. However, while bacterial T-RFLP fingerprinting did not differentiate the east/west dune slopes (**Figure 3**), pyrosequencing analysis detected minor differences (**Figure 4**) which could be artifactual and related to the low number of samples analyzed.

A total of thirty bacterial phyla were observed in the dune samples (Table S1), with only seven (Acidobacteria, Actinobacteria, Bacteroidetes, Chloroflexi, Deinococcus/Thermus, Firmicutes, and Proteobacteria) showing relative abundances

greater than 1% (**Figure 4**). Proteobacteria, Actinobacteria, and Bacteroidetes were the dominant phyla in the 7 dune environments studied (Table S1). Bacterial phyla habitatfiltering could be observed with Proteobacteria dominating sand dune samples (i.e., Top, Slope and Base), ranging from 42 to 49%, and Actinobacteria dominating the Interdune biotope (51%; **Figure 4**). Moreover, the Chloroflexi, Firmicutes, and Acidobacteria (**Figure 4**) as well as the Alpha-, Beta-, Gamma-, and Deltaproteobacteria (**Figure 5**) displayed dune habitat-specific relative abundances: Chloroflexi, Firmicutes, Alpha-, Gamma-, and Delta-proteobacteria presented decreasing abundances from dune Top to the Interdune, while Acidobacteria and Betaproteobacteria presented the opposite trend (**Figure 4**). Similarly, the three most abundant genera of the dominant Actinobacteria (Geodermatophilus sp., Blastococcus sp., and Arthrobacter sp.) and Proteobacteria (Microvirga sp., Massilia sp., and Novosphingobium sp.) phyla displayed dune habitat specific abundances (Table S2). Altogether, these results strongly suggest habitat-filtration to be a determining process defining dune bacterial community assemblies.

This is confirmed by the fact that, apart from the dune Tops, each dune zone presented different and significantly associated (p < 0.05) bacterial indicator species (Dufrêne and Legendre, 1997). Two actinobacterial OTUs were associated with the Interdune: OTU24 from the genus Geodermatophilus sp. and OTU88, a Modestobacter sp. Three actinobacterial genera were also significantly linked to the eastern dune Slope (Kineococcus sp. [OTU73 and OTU77], Arthrobacter sp. [0TU16] and Rathayibacter sp. [0TU46]). Contrastingly, the two west dune slope indicator species were from the Flavisolibacter sp. (OTU6) and Microvirga sp. (OTU23) genera. Finally the dune

Bases (east and west) were characterized by indicator species OTU74 from the Acidobacteria Gp16 class and OTU85 from the genus Kineococcus sp.

# Dune Zone-specific Metabolizing Bacterial Community

FDA hydrolysis was used as a proxy for the active bacterial metabolism (Chrzanowski et al., 1984). Active bacterial communities were detected in all dune areas studied (**Figure 6**), with the highest activity detected in the Interdune samples (average values 0.7–0.75µg fluorescein g−1 soil h−1 ). Lowest activity was observed in Top and Slope samples (0.05–0.15µg fluorescein g−1 soil h−1 ) with Base samples showing intermediate values (0.15–0.4µg fluorescein g−1 soil h−1 ). The results of the Tukey's HSD test on the FDA hydrolysis data indicated that there were only three statistically distinct functional bacterial communities in the dune system: Interdune, Base and Slope/Top (**Figure 6**).

# Identification of Environmental Factors Shaping Dune Bacterial Community Structures

To further investigate relationships between abiotic factors and bacterial community structure (i.e., T-RFLP results), we performed redundancy analysis (RDA; **Figure 7**). RDA indicated that of the 21 soil attributes measured, 15 parameters were significant in explaining bacterial community structure variability (p < 0.01). However, taken together they only explain 21% of this variability. These included mineral [Na, K, P, NH4], metallic [Fe, Mn, Al] and granulometry factors [VCS, CS, MS, FS and (VFS/SC/ST/CY)], and nutrients [P, NH4]. Particle size distribution of the sand dune soils was an important environmental determinant to emerge from the RDA analysis. Our data (**Table 1**) show very clear grain size segregation, with the finest sands, silts and clays predominantly in the Interdune samples and the medium and coarser grained sands on the dunes themselves. The data indicated that fine grains, silts and clays [VFS, SC, ST, CY] are positively correlated with metals [Fe/Mn/Al], minerals [Ca/Mg] and some nutrients (C, P, NO− 3 - N, NH+ 4 -N) (**Figure 7** and Figure S1C).

Two main clusters of bacterial community samples were separated along the RDA axis1 as shown on **Figure 7**: the dune

Base and Interdune bacterial community samples (cluster 1) and the dune Slope and Top samples (cluster 2). The bacterial community structures from cluster 1 were mainly driven by VFS, VCS and NH+ 4 -N, K, Fe, and P contents while those of cluster 2 by MS and Na contents. Communities from both clusters were further separated along the RDA axis2 according to their dune biotope of origin, with granulometry parameters [FS and CS] differentiating the structures of the Interdune and dune Base communities as well as the ones from the dune Top and Slope (**Figure 7**).

Overall, these results demonstrate that environmental gradients in the four defined biotopes (i.e., dune Top, Slope, Base, and Interdune) significantly impact the bacterial community structures of dune/interdune transects.

# Discussion

The geomorphology and the macro-fauna/flora of the 34,000 km2 area of the Namib Sand Sea have been extensively studied (Seely and Pallett, 2008; Livingstone, 2013). However, arguably the most productive component of this depauperate environment, the edaphic microbial community (Makhalanyane et al., 2015), has received little attention at either local (a unique dune system) or regional (multiple dune systems) scales (Jacobson, 1997). Moreover, the drivers of desert dune microbial community assembly generally remain unclear.

This study was designed to evaluate the importance of deterministic processes (e.g., habitat filtering) in influencing bacterial communities in the Namib Sand Sea dune system. We provide a comprehensive analysis of a localized site, with 125 soil samples collected from five parallel 2 km dune/interdune transects. All samples were characterized by 16S rRNA gene T-RFLP fingerprinting and the determination of multiple physicochemical parameters. The use of multiple parallel transects was intended to minimize the effects of local heterogeneity and spatial unevenness, and to generate statistically robust and representative datasets (Prosser, 2010). 16S rRNA gene

amplicons from representative samples were also pyrosequenced to provide more insight into the phylogenetic composition of dune/interdune transects and evaluate if a priori defined niches would present different community assemblages.

# Soil Properties of the Namib Desert Dune/Interdune Transect Define Four Biotopes

The physico-chemical characteristics of soils along the Namib Desert dune-interdune transects clustered into four distinct datasets (**Figure 2**; **Table 1**): dune Top, Slope, Base, and Interdune. Granulometry data were typical for Namib Desert dune systems, with the highest proportion of fine sand on dune summit; a phenomenon commonly attributed to segregation of aeolian deposits from the gravel Interdune zone, and general concentration of coarse sands toward the dune Base (Lancaster, 1981). High Na content in the dune top samples is thought to be the result of salt capture from marine fog, the commonest source of moisture in this particular Namib Desert environment (Abrams et al., 1997; Hachfeld, 2000; Eckardt et al., 2013b). The high metal content recorded in the Interdune zone is likely to reflect the underlying geology, which includes red Tsondab sandstone (Besler, 1996), but may also be a consequence of the higher proportion of the fine soil particle fractions (VFS, SC, ST, and CY), the presence of which reduce percolation and enhance mineral accumulation and persistence (Chamizo et al., 2012).

# Bacterial Community Structures, Diversities and Potential Function Suggest Habitat Filtration in Dune/Interdune Transects

A total of 30 bacterial phyla were detected in Namib Desert dune soils and all samples were dominated by Proteobacteria, Actinobacteria, and Bacteroidetes, which has also been observed in other desert systems (Makhalanyane et al., 2015). Namib Desert dune samples showed similar phylogenetic profiles to dune systems from the Gobi (Mongolia) and Taklamaken (China) deserts (An et al., 2013), but with lower OTU diversity than dune slope and interdune of the Gurbantunggut (China) desert (Li et al., 2015). The latter, however, was characterized by cyanobacteria-rich biological soil crusts (BSC), which were not observed in the Namib Sand Sea at the time of sampling.

Both 16S rRNA gene T-RFLP (**Figures 3**, **7**; Figure S1A) and pyrosequencing (**Figures 4**, **5**; **Table 1** and Table S1; Figure S1B) data strongly support the concept of niche-partitioning/habitatfiltering in the Namib Desert dune/interdune system (Dumbrell et al., 2010), as both the structures and phylogenetic compositions of the bacterial communities demonstrated similar clustering than those of the physico-chemical characteristics (**Figure 2**). This was further supported by the low number of OTUs (15% from T-RFLP and 2% from pyrosequencing) shared by all dune biotopes (Figure S1). Taken together, these results strongly imply that habitat filtering, possibly via environmental physico-chemical gradients, is a major driver of bacterial community composition in dune ecosystems (Brankatschk et al., 2013; Liu et al., 2014).

While bacterial community T-RFLP fingerprinting and pyrosequencing showed highly comparable overall microbial community patterns (**Figures 3**, **4**), pyrosequencing data showed minor differences in bacterial community structures between the east and west dune slopes (**Figure 4**). As the edaphic characteristics of these dune aspects were not significantly different, it is suggested that other environmental parameters, such as aeolian disturbance (Thomas and Dougill, 2007) or solar radiation (Jacobs and Sundin, 2001) differentials, could account for these differences.

# Dune/Interdune Transect Habitat-filters

Distance based RDA analysis showed that 15 edaphic parameters, including mineral [Na, K], metallic [Fe/Mn/Al] and granulometry factors [VCS, CS, MS, FS and (VFS/SC/ST/CY)], and nutrients [P, NH4] were strongly correlated with bacterial community structure variability in dune soil environments (**Figure 7**, **Table 1**). Using Na as a proxy for soil salinity, we conclude that this factor plays an important role in directing dune microbial community compositions. Variations in salinity have previously been shown to be a strong bacterial habitat filter in Namib Desert gravel plain soils and hypoliths (Stomeo et al., 2013). The importance of soluble salts in dictating bacterial species abundance is supported by the observation that halotolerant microorganisms are readily isolated from interdune soils in the Negev Desert (Yu and Steinberger, 2012b) and Thar Desert sands (Sharma et al., 2013).

The role of K as a habitat filter may be related to intracellular osmoregulation (Fierer et al., 2012) while that of metals [Fe/Mn/Al] could be linked to their bioavailability, typically controlled by solubility and/or adsorption onto mineral surfaces (Cornell and Schwertmann, 2003). In this resource-limited desert environment, it was unsurprising to identify the nutrients P and NH4 as environmental drivers of dune bacterial community structures. The bioavailability of phosphorus and nitrogen has indeed previously been found to impact desert soil microbial functional and structural properties (Austin et al., 2004; Bell et al., 2014). We note that the use of P-Bray I to quantify total P underestimated the P contents in these mildly alkaline soils (Bowman and Vigil, 2002). In future studies, a more appropriate method such as the sodium bicarbonate procedure (Olsen et al., 1954) should be used to properly evaluate the bioavailable P contents in these desert soils and its potential as a factor shaping dune bacterial communities. The precise role of these dbRDAdetermined environmental parameters in shaping dune bacterial community structures and functions could be further tested in controlled experiments (e.g., micro- or mesocosms; Lan et al., 2014).

The single most important environmental determinant to emerge from these data is the desert dune soil particle size distribution. Specific microbe-particle associations have already been demonstrated (e.g., Zhang et al., 2007; Carson et al., 2010). This could be related to the fact that differences in soil structure can manipulate soil nutrient status (Zhang et al., 2007; Carson et al., 2010). Our data do indicate that fine grains, silts and clays [VFS, SC, ST, CY] are positively correlated with metals [Fe/Mn/Al], minerals [Ca/Mg] and some nutrients (C, P, NO− 3 -N, NH+ 4 -N) (**Figure 2**). This suggests that the relationship between bacterial diversity, soil structure and nutrient contents needs further investigation in dune soils as subtle differences in sand mineralogy and micronutrient contents has been shown to significantly impact edaphic community structuring (Welz et al., 2014).

We noted that 79% of the variation in soil bacterial community structures remains unexplained (**Figure 5**), although this is not uncommon in microbial ecology studies (e.g., Cottenie, 2005; Dumbrell et al., 2010). The high proportion of unaccounted variation might be attributed to other deterministic environmental factors not recorded in this study (e.g., wind transport, effect of sand burial stress or changes in humidity, and/or other soil physicochemistries) or the effects of microbial interactions (Vellend, 2010; Caruso et al., 2011). We also cannot rule out the role of stochastic processes in the assembly of dune communities, as, at a global scale, it has been shown significantly drive bacterial community assembly in deserts (Caruso et al., 2011).

The Namib Desert sand dunes are physically dynamic environments with high levels of wind erosion and dispersion (Eckardt et al., 2013a; Livingstone, 2013). It could be argued that in apparently well-connected habitats such as those of the Namib Sand Sea, dispersal rates should be sufficiently high so as to lead to homogenization of microbial communities due to mass-mixing effects (Leibold and Norberg, 2004). However, our bacterial fingerprint results and phylogenetic analyses show that only a limited number of OTUs were shared across all four dune habitats and that the number of OTUs unique to each area and OTUs shared by contiguous habitats was also low (Figure S1). This suggests that aeolian transport and cellular dispersal is not a substantial driver of sand dune bacterial community structure. However, it cannot be completely excluded. Indeed, in the dune Top zone, which experiences the highest wind disturbance (Eckardt et al., 2013a), the lowest levels of bacterial activity were observed (as determined by fluorescein diacetate hydrolysis: **Figure 6**) and no indicator species was identified. It could be argued that this zone does not represent a stable community and that the detected activity could be derived from microorganisms bound to wind-transported mineral grains (Yamaguchi et al., 2012).

# Phylogenetic Composition of Namib Desert Dune Biotope Soils in Relation to Pedogenesis and Ecosystem Development

Many studies have shown that microbial communities evolve with the development of ecosystems, over short (weeks) or long (years/hundreds of years) temporal scales (e.g., Ramond et al., 2012; Ferrenberg et al., 2013; Jangid et al., 2013). For example, in a 6500 year dune soil chronosequence, Betaproteobacteria, Bacteroidetes, Actinobacteria, Cyanobacteria, and Firmicutes were shown to dominate "young soil ecosystems"; i.e., recent or often-disturbed soil environments (Jangid et al., 2013). Contrastingly, the relative abundances of Acidobacteria, Alphaproteobacteria, and Planctomycetes were higher in "aged soils"; i.e., more structured and stabilized environments (Jangid et al., 2013). Similarly, continuous physical disturbance has been shown to influence community composition in cold desert biocrust soils (Kuske et al., 2012), with continuously trampled soils showing decreasing cyanobacterial abundances with corresponding increases in Actinobacteria, Chloroflexi, and Bacteroidetes, compared to undisturbed soils.

As previously observed in a 6500 year dune soil chronosequence (Jangid et al., 2013), Betaproteobacteria, Bacteroidetes, and Firmicutes were more abundant in the more dynamic dune Tops and Slopes, while the relative abundances of Acidobacteria and Alphaproteobacteria were higher in the more stable Interdune soils. Cyanobacteria were almost exclusively detected in the stable Interdune soils (albeit in proportionally very low percentages: 0.46%), which were dominated by Actinobacteria, comprising over 50% of the interdune community (Table S1). Our results are remarkably consistent with the comparative data obtained from disturbed and stable biocrust systems (Kuske et al., 2012) and suggest that the degree of stability or perturbation of the Namib Desert soils may be one of the key deterministic factors contributing to the high proportion of unaccounted variation in community structure of the Namib Desert dune systems.

# Concluding Remarks

Our results clearly showed that, in dune/interdune transects of the Namib Desert, soil physico-chemistry define four dune zones, namely the dune Tops, Slopes, Bases and Interdune. Furthermore, we demonstrated that each of these zones presented unique bacterial communities. Altogether, this strongly suggested that these communities were selected through habitat filters, i.e., environmental deterministic factors (Dumbrell et al., 2010). The measurement of 21 environmental factors enabled us to define that habitat stability, soil texture and mineral and nutrient contents were influencing dune/interdune transect bacterial communities. However, almost 80% of the dune/interdune bacterial community's variability remained unexplained indicating that other unmeasured deterministic factors or processes, such as species' interaction, dispersal and/or stochasticity (Vellend, 2010), may be involved. Consequently, to further investigate bacterial community assembly and its driver(s) in dune systems, such study should be extended to multiple dune/interdune environments.

# Author Contributions

SR, JR, and DC designed the experiment. SR performed the lab work and analyzed the data with some help from JR. JR, BJ, SR, MS, and DC participated in writing the manuscript.

# References

DC contributed to the reagents, materials and analysis tools required. MS' field knowledge was necessary in the design experiment.

# Acknowledgments

The authors thank the University of Pretoria Genomics Research Institute and the National Research Foundation (NRF) of South Africa for financial support. We also thank all the members of CMEG who assisted with fieldwork in the Namib Desert. Finally, we thank the reviewers as for their very constructive remarks that improved the manuscript.

# Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fmicb. 2015.00845

**Conflict of Interest Statement:** 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.

Copyright © 2015 Ronca, Ramond, Jones, Seely and Cowan. This is an openaccess article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Soil bacterial and fungal community responses to nitrogen addition across soil depth and microhabitat in an arid shrubland

#### Rebecca C. Mueller 1 , Jayne Belnap2 and Cheryl R. Kuske1 \*

1 Bioscience Division, Los Alamos National Laboratory, Los Alamos, NM, USA, 2 Southwest Biological Science Center, United States Geological Survey, Moab, UT, USA

#### Edited by:

Thulani Peter Makhalanyane, University of Pretoria, South Africa

#### Reviewed by:

Ulas Karaoz, Lawrence Berkeley National Laboratory, USA Kimberley Warren-Rhodes, SETI Institute and NASA Ames Research Center, USA Burkhard Buedel, University of Kaiserslautern, Germany

#### \*Correspondence:

Cheryl R. Kuske, Bioscience Division, Los Alamos National Laboratory, M888, Los Alamos, NM 87545, USA<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Terrestrial Microbiology, a section of the journal Frontiers in Microbiology

Received: 23 June 2015 Accepted: 17 August 2015 Published: 04 September 2015

#### Citation:

Mueller RC, Belnap J and Kuske CR (2015) Soil bacterial and fungal community responses to nitrogen addition across soil depth and microhabitat in an arid shrubland. Front. Microbiol. 6:891. doi: 10.3389/fmicb.2015.00891 Arid shrublands are stressful environments, typified by alkaline soils low in organic matter, with biologically-limiting extremes in water availability, temperature, and UV radiation. The widely-spaced plants and interspace biological soil crusts in these regions provide soil nutrients in a localized fashion, creating a mosaic pattern of plant- or crust-associated microhabitats with distinct nutrient composition. With sporadic and limited rainfall, nutrients are primarily retained in the shallow surface soil, patterning biological activity. We examined soil bacterial and fungal community responses to simulated nitrogen (N) deposition in an arid Larrea tridentata-Ambrosia dumosa field experiment in southern Nevada, USA, using high-throughput sequencing of ribosomal RNA genes. To examine potential interactions among the N application, microhabitat and soil depth, we sampled soils associated with shrub canopies and interspace biological crusts at two soil depths (0–0.5 or 0–10 cm) across the N-amendment gradient (0, 7, and 15 kg ha−1 yr−1 ). We hypothesized that localized compositional differences in soil microbiota would constrain the impacts of N addition to a microhabitat distribution that would reflect highly localized geochemical conditions and microbial community composition. The richness and community composition of both bacterial and fungal communities differed significantly by microhabitat and with soil depth in each microhabitat. Only bacterial communities exhibited significant responses to the N addition. Community composition correlated with microhabitat and depth differences in soil geochemical features. Given the distinct roles of soil bacteria and fungi in major nutrient cycles, the resilience of fungi and sensitivity of bacteria to N amendments suggests that increased N input predicted for many arid ecosystems could shift nutrient cycling toward pathways driven primarily by fungal communities.

Keywords: dryland, shrubland, soil fungal community, soil bacterial community, global change, biological soil crusts, microhabitat, ribosomal RNA

# Introduction

Nitrogen (N) deposition due to human activities has increased dramatically and is predicted to rise with human population growth (Galloway et al., 2008). The effects of N additions on soil microbial communities have been well-studied in mesic forests in the eastern U.S., but the effects of adding N to dryland ecosystems are poorly understood. Arid lands have distinct characteristics compared to forests, with multiple limiting nutrients, including water, N, and P, and commonly have low inputs of carbon (Collins et al., 2008; Sterflinger et al., 2012), and it is unlikely that microbial responses to N can be accurately inferred from patterns documented in more mesic systems, particularly given documented interactions between soil properties and N responses. Historically, arid regions have received lower natural N inputs (∼1–2 kg ha−1 yr−1 ; Phoenix et al., 2006), but it is not known if these lower levels have similar impacts comparable to higher inputs to more mesic regions. In addition, the native microbial communities in drylands are highly dissimilar from those found in mesic environments (Fierer et al., 2003; Dunbar et al., 2012). The effects of experimental climate change on soil microbial communities can depend upon the composition of the native community (Hawkes and Keitt, 2015), suggesting that variation among ecosystems could alter patterns of community structure and response to N deposition. Urban areas are sources of atmospheric N and are increasing in arid ecosystems (Galloway et al., 2004); thus, a future increase in N deposition is likely in these regions. Understanding soil microbial community responses to additional N input will provide insight on how we may expect the functions of dryland soils to change in response to future conditions, and how we may effectively manage nitrogen inputs into drylands to maintain productivity and resilience of microbial communities.

Drylands are typified by high heterogeneity in seasonal climate, the distribution of nutrients across the landscape, and physical microhabitats. Many hot dryland ecosystems experience strong seasonal shifts, with pronounced dry and rainy seasons and large fluctuations in temperature (Collins et al., 2008). These systems are also characterized by widely spaced vegetation, with biological soil crusts (biocrust) occupying the soil surface between the vascular plants (Belnap, 2003). Habitat heterogeneity has been shown to promote beta diversity across small geographic ranges for plants, insects, and microorganisms (Amarasekare, 2003), where differences in climatic factors are minimal. Previous studies of microbial community biomass and composition in dryland systems suggests these features are strongly correlated with differences in microhabitat (Garcia-Pichel et al., 2003; Schade and Hobbie, 2005; Bates and Garcia-Pichel, 2009; Steven et al., 2013, 2014), which can, in turn, shift function and subsequently alter microbial responses to changes in climate such as altered precipitation (Johnson et al., 2012; Delgado-Baquerizo et al., 2014). In addition, at a broad scale, microbial metabolic pathways have also been shown to differ between rhizosphere and biocrust soils (Steven et al., 2014).

Although most dryland soils do not have highly weathered soils with characteristic organic horizons, soil nutrients are highly stratified by depth (Pointing and Belnap, 2012). In addition, drylands are characterized by highly patchy distribution of vegetation (Collins et al., 2008), leading to high levels of heterogeneity at both vertical and horizontal scales. A full understanding of the microbial response and ecosystem effects of N addition in arid systems therefore requires examination of multiple microhabitats and soil depths. Strong vertical stratification of communities across shallow soil horizons has also been found for bacterial communities in dryland ecosystems (Kuske et al., 2002; Steven et al., 2013), with significant differences in communities between plant and crust-associated soils (Steven et al., 2013), but next generation sequencing technologies now provide the means for in-depth analysis of microbial community shifts. Using high throughput sequencing, we examined fungal and bacterial community response to N addition within two microhabitats: the plant root zone and interspaces soils covered by biocrusts. Soils were sampled at two depths in each microhabitat (0–0.5 cm surface veneer, and 0–10 cm soil).

Previous studies have documented patterns of bacterial or fungal composition in arid lands (e.g., Kuske et al., 2002; Bates et al., 2012; Steven et al., 2013, 2014), but have not examined the effects of N. In addition, few studies have utilized high-throughput sequencing methods to characterize microbial communities within drylands (Makhalanyane et al., 2015). Correlated shifts in bacterial communities, catabolic profiles, and metagenomic surveys have been shown in experimental N sites (Amarasekare, 2003; Fierer et al., 2012), but with no concurrent analysis of the fungal community. Although bacterial biomass generally exceeds that of fungal biomass in arid ecosystems (Bates and Garcia-Pichel, 2009), the majority of biogeochemical cycles in terrestrial ecosystems are driven by both bacteria and fungi (van der Heijden et al., 2008), so an understanding of the responses of both bacterial and fungal communities is needed to predict long-term ecosystem effects of N inputs. Consequently, we conducted compositional analysis of these two communities, combined with examination of nutrient transformations (Sinsabaugh et al., 2015), to provide insights into both the community- and ecosystem-level effects of N deposition on drylands.

# Materials and Methods

The study site is located within the Lake Mead National Recreation Area in southern Nevada (36.009 N, 114.797 W). Two species of shrubs, Larrea tridentate and Ambrosia dumosa comprise the dominant vegetation, interspersed with patches of poorly developed biological soil crusts (biocrusts). Soils are alluvial, with high sand content (80%), low organic carbon (0.35%) and slightly basic pH (7.4 in bulk soils). Mean high and low annual temperatures are 23.8◦C and 11.1◦C, respectively, with mean annual precipitation of 19.8 cm based on data collected from a nearby weather station.

Within a 100 × 100 m plot, fifteen 2 × 2 m subplots were established, centered around similarly-sized Ambrosia dumosa shrubs. To examine the effects of N on plant and soil microbial communities, N was experimentally manipulated at two levels: 7 and 15 kg−1 ha−1 yr was added, along with an ambient control, hereafter referred to as N7, N15 and ambient. Nitrogen was applied in the form of ammonium nitrate in water every 2 months from March 2012 to May 2013, when soils were harvested.

To examine potential impacts of N amendment on bacterial and fungal community composition, we sampled soils from the two microhabitats (shrub and interspace) across two soil depths (surface and sub-surface) (**Figure 1**). Surface soils were

sampled to a depth of 0.5 cm, and sub-surface soils were sampled to a depth of 10 cm. Five replicates per N treatment were collected from four sampling locations: shrub canopy surface at 0–0.5 cm depth (canopy biocrust), shrub canopy bulk soils at 0–10 cm depth (canopy bulk), interspace surface at 0–0.5 cm depth (interspace biocrust), and interspace bulk soil at 0–10 cm depth (interspace bulk). A total of 60 samples were collected. Soils were homogenized by hand and total nucleic acids were extracted from 0.5 g of soil using the MP Biomedical FastDNA for Soil kit (MP Biomedical, Santa Ana CA, USA) according to the manufacturer's instructions.

Bacterial and fungal communities were targeted for Illumina high-throughput sequencing using a two-stage PCR approach described previously (Mueller et al., 2014). The V3-V4 region of the bacterial 16S rRNA gene was amplified using the F515- R806 primer pair (Bates et al., 2010), and the D2 hyper-variable region of the fungal LSU rRNA gene was amplified using the newly designed LR22R primer (TAGCGMACAAGTASMGTGA) with the LR3 primer as the reverse (http://sites.biology.duke.edu/ fungi/mycolab/primers.htm). Samples were barcoded to facilitate multiplexing using a combinatorial approach (Gloor et al., 2010), with a unique 6 bp barcode inserted into both the forward and reverse primer. Illumina-specific sequences were added in a second PCR. Amplicons were cleaned using the MoBio –htp UltraClean 96-well plate clean up kit according to the manufacturer's directions, but with a modification to reduce the amount of binding buffer from 5X to 3X sample volume to facilitate removal of longer primer-dimers. Individual samples were quantified using the Quant-It high-sensitivity double stranded DNA kit (Invitrogen, Carlsbad CA USA) using the BioTek Synergy H1 microplate reader and combined at equimolar concentrations. The multiplexed library was cleaned using the MoBio UltraClean PCR clean up kit and eluted in 50µl EB buffer. Dissimilar values from fluorometry and Illumina-specific qPCR indicated that the second PCR step did not add the Illumina-specific primers to all samples, and a final PCR, with conditions identical to the second PCR, was performed on the library, for a total of three cycles. Library length and concentration was verified using a bioanalyzer and qPCR, respectively. Amplicons were sequenced using paired end 250 bp reads on the Illumina MiSeq platform at Los Alamos National Laboratory. To offset the low diversity of nucleotides in the initial sequence common in amplicon sequences, genomic DNA from multiple bacterial isolates was added at approximately 30%.

Bacterial and fungal community sequence datasets were de-multiplexed individually using QIIME (Caporaso et al., 2010b), with quality filtering to remove any sequence with a mismatch to the barcode or primer sequence. Downstream sequence processing was conducted using UPARSE (Edgar, 2013). Sequences within an expected error greater than 0.5 and singleton sequences were removed. Putative chimeras were identified de novo using UCHIME (Edgar et al., 2011) and removed, and OTUs were delineated at 97% sequence similarity. The OTU table was exported, and all statistical analyses were conducted using the R statistical platform (r-sourceforge.org). Representative fungal and bacterial sequences were selected and classified against their respective Ribosomal Database Project (RDP) databases using the online tool. Any sequence identified as non-Fungal or non-Bacterial was removed from subsequent analyses. Sequences were deposited to the MG-RAST server (Meyer et al., 2008) under IDs 4635375.3 (Bacteria) and 4635376.3 (Fungi).

To construct phylogenetic trees, we utilized reference phylogenies constructed using partial to full-length ribosomal sequences. For the bacteria, we used the tree from Kembel et al. (2012), and for fungi we used a tree constructed from fungal sequences generated by the Assembling the Fungal Tree of Life Project (Celio et al., 2007). Representative OTU sequences were aligned to the reference taxa using PyNAST (Caporaso et al., 2010a) and placed on their respective reference phylogenies using pplacer (Matsen et al., 2010).

Statistical analyses were conducted using the R libraries picante (Kembel et al., 2010), vegan (Oksanen et al., 2013), and phyloseq (McMurdie and Holmes, 2013). Analyses were conducted by randomly selecting 5200 sequences for bacteria and 4500 sequences for fungi, with 99 trials for each analysis to limit bias resulting from random sampling of communities. Three bacterial and two fungal samples had insufficient sequencing depth and were not included in the analyses. Richness measures, including OTU richness and phylogenetic diversity (PD) (Faith, 1992) were compared using Three-Way ANOVA with N treatment, microhabitat, and soil depth as fixed factors. Community composition based on OTU similarity (Bray-Curtis dissimilarity metric) and phylogenetic similarity using the weighted UniFrac distance (Lozupone and Knight, 2005) were compared using PERMANOVA with all three factors included in the model (Anderson, 2001). In addition, to examine responses at broader taxonomic scales, we grouped fungal and bacterial OTUs at class and phylum levels, respectively, and tested for significant differences in relative abundance across N treatment, soil depth, and microhabitat using Three-Way ANOVA.

# Results

Following removal of low-quality sequences, a total of 1,340,717 fungal LSU sequences and 1,466,336 bacterial 16S sequences were included in the analysis. Across all samples, we observed 1990 fungal OTUs and 10167 bacterial OTUs. Bootstrap estimates of gamma diversity indicated that both communities were slightly under-sampled; the bootstrap estimate was 2226 (±45 SE) for fungi and 10884 (±160 SE) for bacteria (**Figure 2**).

# Effects of Microhabitat and Soil Depth

Microhabitat and soil depth each had strong effects on fungal and bacterial community richness and PD. For the fungal community, soil depth was the strongest driver of diversity

differences for both OTU richness (F = 128, p < 0.001) and PD (F = 11.7, p = 0.001). Surprisingly, richness and PD were higher in the bulk soil depth (**Figure 3**). OTU richness did not vary significantly between the shrub canopy and the interspace (F = 3.34, p = 0.07), although there was a trend for higher fungal richness in association with interspace compared to shrub canopy soils (**Figure 3**). However, PD of fungal communities associated with the interspaces was significantly higher than those found in the shrub canopy (F = 4.42, p = 0.04).

The bacterial community also showed strong differences in relation to soil depth. As with the fungal community, OTU richness (F = 140, p < 0.001) and PD (F = 137, p < 0.001) were significantly higher in the deeper bulk soils (**Figure 3**). Bacterial richness also differed significantly between the two microhabitats, with higher richness (F = 0.24, p < 0.001) and PD (F = 58, p < 0.001) observed in the interspaces than in the shrub canopy when analyzed across both soil depths.

Although for fungal and bacterial communities, the strongest richness differences were found between the two soil depths, the largest driver of variation in community composition was microhabitat (**Figures 4**, **5**). This pattern was observed using both OTU-based measures (Fungi: F = 8.89, R2 = 0.13, p = 0.001, Bacteria: F = 6.10, R2 = 0.10, p = 0.001) and phylogeneticbased measures (Fungi: F = 11.6, R2 = 0.15, p = 0.001, Bacteria: F = 23.9, R2 = 0.23, p = 0.001). Both communities were also delineated by soil depth, but this factor explained a lower proportion of community variation for both OTU-based (Fungi: F = 2.50, R2 = 0.04, p = 0.001, Bacteria: F = 2.62, R2 = 0.04, p = 0.03) (**Figure 4**) and phylogenetic measures (Fungi: F = 9.49, R2 = 0.12, p = 0.001, Bacteria: F = 15.7, R2 = 0.15, p = 0.001) of community similarity (**Figure 5**).

### Nitrogen Effects on Microbial Communities

We found no evidence for fungal community responses to N addition. We observed no shifts in OTU richness in response to N addition (F = 0.48, p = 0.62) and no change in Faith's PD (F = 0.11, p = 0.90) (**Figures 4**, **5**). Neither OTU-based community composition (F = 1.17, R2 = 0.04, p = 0.18), nor phylogenetic composition (F = 1.18, R2 = 0.04, p = 0.27) shifted in response to experimental N additions (**Figures 4**, **5**).

In contrast to responses in the fungal community, we found consistent negative responses to N in the bacterial community diversity indices (**Figure 3**). Both richness (F = 9.85, p < 0.001) and PD (F = 11.7, p < 0.001) declined in N-amended plots. Tukey's Honestly Significant test indicated that for both measures, intermediate and high levels of N (N7 and N15) were not significantly different (p > 0.2 for each). Nitrogen addition was associated with community shifts using both OTUbased (F = 1.34, R2 = 0.04, p = 0.02) and phylogenetic measures (F = 3.42, R2 = 0.07, p = 0.001; **Figure 4**). Using the UniFrac distance measure, we also found a significant N X soil horizon depth interaction (F = 1.87, R2 = 0.04, p = 0.04), where community shifts were more pronounced at lower soil depths.

## Broad Taxonomic Responses

The phylum Ascomycota was the dominant fungal phylum within this shrubland. At finer taxonomic scales, the most

abundant fungal classes included the Dothideomycetes, Pezizomycetes, and the Basidiomycete class Agaricomycetes. The bacterial community was dominated by the Proteobacteria, followed by the Actinobacteria and Acidobacteria (**Figure 6**).

When bacterial and fungal communities were examined by grouping OTUs into bacterial phyla or fungal classes, respectively, we found that only a single fungal class, the Sordariomycetes, showed a response to N addition, with a nearly linear increase in relative abundance in N amended plots (F = 4.34, p = 0.02). Three classes, the Leotiomycetes, Sordariomycetes, and Eurotiomycetes, were all significantly more abundant at deeper soil depths (for each, p < 0.03). No significant effect of microhabitat on the abundance of any fungal class was observed.

In contrast, numerous bacterial phyla showed significant shifts in relative abundance in response to N additions (**Table 1**). The two most abundant groups, the Actinobacteria and Proteobacteria, showed negative responses at the highest level of N input. Responses for the classes within the Proteobacteria were consistent, with declines observed for alpha, beta, delta, gamma, and epsilon classes. Four other bacterial phyla, the Acidobacteria, Chloroflexi, Firmicutes, and Verrucomicrobia, showed significant positive responses to N amendment. In addition, a number of bacterial phyla showed significant differences across soil depth and microhabitat type. The relative abundance of most phyla was higher at lower soil depths, with the exception of Proteobacteria and Verrucomicrobia, which exhibited higher relative abundance in the surface soils. Microhabitat also affected a large number of bacterial phyla. For example, the relative abundance of Acidobacteria was higher in the interspace soils, while the proportion of Bacteroidetes was higher in shrub-associated soils (**Table 1**).

FIGURE 4 | Bray-Curtis dissimilarity among fungal or bacterial communities, displayed as 2D MDS plots. Plots are color coded to show microhabitat (red), soil depth (green), or N amendment (blue) differences. While results are presented for individual treatments, results were analyzed using a three-way framework with PERMANOVA. See Results text for statistical support.

# Discussion

Our study objective was to determine the cumulative impacts of N deposition, delivered every 2 months for 14 months, on the microbial communities within an arid shrubland ecosystem. Because of its mosaic nature, we partitioned the landscape into compartments and separately assessed the N amendment impacts on those compartments. This study represents the first high-depth-of-coverage, taxonomic survey of both the fungal and bacterial communities in four shrubland microhabitats—shrub canopy and interspace soils, each stratified at two shallow depths. We also examined the sensitivity of those fungal and bacterial communities to short-term (1 yr), relatively low concentrations of N amendment.

# Fungal and Bacterial Community Differences with Soil Depth

In a companion report to this study, bacterial and fungal biomass, measured by qPCR, was significantly higher in the 0–0.5 cm biocrusts than in 0–10 cm depth bulk soil (>100% and >200% higher, respectively; Sinsabaugh et al., 2015). This is a typical pattern for arid shrublands and grasslands, where, by multiple

measures, microbial biomass has been found to be concentrated in the upper few cm of soil (Dunbar et al., 2012; Pointing and Belnap, 2012; Steven et al., 2013). Biocrust biomass may differ substantially with soil type (Steven et al., 2012) and with biocrust developmental stage (Garcia-Pichel et al., 2003; Yeager et al., 2004), particularly in Cyanobacterial biocrusts where the proportion and density of the dominant cyanobacteria in the biocrust controls the overall biomass. Our study partitioned soil depth into intervals at the cm scale; however, biomass of biocrust soils is stratified even at the mm interval (Garcia-Pichel et al., 2003).