Soil heterogeneity and the distribution of native grasses in California: Can soil properties inform restoration plans?

Abstract

When historical vegetation patterns are unknown and local environments are highly degraded, the relationship between plant species distributions and environmental properties may provide a means to determine which species are suitable for individual restoration sites. We investigated the role of edaphic variation in explaining the distributions of three native bunchgrass species (Bromus carinatus, Elymus glaucus and Nassella pulchra) among central California mainland and island grasslands. The relative contribution of soil properties and spatial variation to native grass species abundance was estimated using canonical redundancy analysis, with subsequent testing of individual variables identified in ordination. Soil variables predicted a significant proportion (22–27%) of the variation in species distributions. Abiotic soil properties that drive species distributions included serpentine substrates and soil texture. Biotic properties that correlated with species distributions were ammonium and nitrogen mineralization rates. Spatial autocorrelation also contributed to species presence or absence at each site, and the significant negative autocorrelation suggested that species interactions and niche differentiation may play a role in species distributions in central California mainland and island grasslands. We explored the application of plant‐environment relationships to ecological restoration for species selection at locations where degradation levels are high and historical communities are unclear.

Introduction

In the face of widespread habitat destruction and biological invasion, land management often requires active restoration of native plant communities. Reintroduction of target species is common practice with the goal of restoring ecosystem composition and, consequently, structure and function (Falk et al. 1996, Chapin et al. 1997). However, reestablishing suites of species—and even single species—can be challenging if site conditions have been greatly altered. Many revegetation efforts fail, and restoration practitioners recognize that success rates can be improved with remediation of physical and biotic environments (Fahselt 2007, Drayton and Primack 2012). Nevertheless, landscape remediation is often inadequate because reference sites are either difficult to identify or nonexistent (Halle and Fattorini 2004), and many sources of disturbance, such as exotic species, represent permanent changes to ecosystem properties (Vitousek 1990).

If we aim to increase the likelihood of plant establishment at restoration sites, we need a better understanding of the factors that drive species distributions in existing environments (Hobbs and Harris 2001). Plant‐environment associations are well‐documented, and while species distributions are primarily controlled by climate at the continental scale, edaphic factors frequently determine species and plant communities at regional or local scales (Cramer and Leemans 1992, Prentice et al. 1992). Examples include nutrient status, pH, and salinity as well as the presence of unusual substrates such as mine tailings and serpentine soils (e.g., Kruckeberg 1951, Gates et al. 1956, Goldberg 1982, Pregitzer and Barnes 1982). Thus, soil properties may serve to predict which species are suitable for reintroduction in cases where historical reference vegetation is unknown, and local environments are highly degraded (Allen and Wilson 1991).

Species distributions reflect spatial processes as well as adaptation to heterogeneous environments (Kramer et al. 2011). Species with limited dispersal show strong relationships with geographic distance, and their distribution declines rapidly with declining environmental similarity. Species with wide dispersal have stronger relationships with ecological distance, and will persist at sites with similar plant communities and environmental conditions at larger geographic scales (Soininen et al. 2007). Much of the current research investigates the contribution of spatial and environmental variation to species range limits and niche availability (e.g., Pan et al. 1998, McCrea et al. 2001, Zhang et al. 2011). Fewer studies have applied plant‐environment and spatial relationships to ecological restoration with the goal of predicting suitable sites for revegetation using selected plant taxa (Allen and Wilson 1991, Corry and McEachern 2000, Volis et al. 2011). These studies offer the potential to improve restoration outcomes via a better understanding of the factors that drive species' distributions.

Dispersal limitations and landscape heterogeneity in soil properties may be important considerations for restoring California grasslands (Harrison et al. 2003, Grman and Suding 2010, Seabloom 2011). Since the period of European colonization, introductions of exotic species have resulted in the invasion of annual Mediterranean grasses and forbs across the state (Bartolome et al. 1986). Native perennial bunchgrass species remain relatively widespread, but historical grassland communities no longer exist and their original composition is unclear (Hamilton 1997, Holstein 2001). Weed control and restoration efforts are underway, but non‐native plants dominate contemporary grasslands, and evidence suggests that ongoing anthropogenic disturbance promotes their persistence (HilleRisLambers et al. 2010, Vallano et al. 2012).

We investigated the role of edaphic variation and spatial processes in explaining the distributions of native bunchgrass species among central California mainland and island grasslands. Our study focused on three perennial bunchgrasses: Bromus carinatus Hook & Arn. (California brome), Elymus glaucus Buckley (blue wild rye), and Nassella pulchra (A. Hitchc.) Barkworth (purple needlegrass) (Hickman 1993). These species are the targets of numerous restoration programs and evidence suggests that they were historically common in coastal woodlands and grasslands (Holstein 2001, Bartolome et al. 2004, Rein et al. 2007). The contribution of spatial and environmental components to bunchgrass species abundance was estimated using canonical analysis and variation partitioning, with subsequent testing of individual variables identified in ordination. Our objectives were to: (1) test the relative importance of spatial structure and soil variables for species distributions, and (2) determine whether soil biotic and abiotic properties predict the potential for target species to reestablish at restoration sites.

Methods

Site characteristics and study species

We studied the distributions of three native perennial bunchgrass species relative to soil factors among two central California mainland sites and five Channel Islands (Fig. 1). Mainland sites included the 2358‐ha University of California Sedgwick Reserve in Santa Ynez, California and Vandenberg Air Force Base near Lompoc, California. Island sites included four northern Channel Islands (Anacapa, Santa Cruz, Santa Rosa and San Miguel) in Channel Islands National Park, and one southern Channel Island (Santa Catalina) where we collected samples within the Catalina Island Conservancy. Soils across the region have a diverse volcanic and sedimentary origin (Norris 2003) and serpentine rock outcrops occur at mainland sites and Santa Catalina Island. Serpentine soils are low‐nutrient substrates, and have a high ratio of magnesium to calcium as well as elevated levels of nickel and chromium (Whittaker 1954, Schechter and Bruns 2008). The regional climate is Mediterranean with hot, dry summers and cool, wet winters (Schoenherr et al. 1999).

Distribution of sampling points among the two California mainland locations, VAFB and Sedgwick Reserve, and the five Channel Islands.

Following more than two centuries of exotic plant and animal introductions, the Channel Islands have been the subject of intense efforts to remove non‐native herbivores such as sheep, pigs, goats and cattle (Donlan et al. 2003). However, introduced plant species are widespread and continue to present challenges for island restoration (Halvorson 1994). On the mainland, livestock have been removed from Sedgwick Reserve, but are still present at VAFB. Mainland grasslands are also targets for restoration and are heavily impacted by exotic plant species. Native bunchgrass communities throughout the region grow in a matrix of introduced, European annual grasses that include Avena spp., B. diandrus, B. hordeaceous, Hordeum murinum, Lolium multiflorum, and Vulpia myuros (Jackson 1985). Other common, introduced species include yellow starthistle (Centaurea solstitialis), fennel (Foeniculum vulgare) and iceplant (e.g., Carpobrotus edulis) (D'Antonio 1993, Bell et al. 2008, Knapp et al. 2009).

Bromus carinatus, E. glaucus, and N. pulchra are wind‐pollinated, non‐rhizomatous bunchgrasses. All three species have awned seeds that are dispersed passively by wind, with some potential for dispersal by animals. Investigation of the factors influencing the persistence of native species in California grasslands has indicated these three bunchgrass species are propagule limited, and seedling success is affected by spatial heterogeneity (Seabloom et al. 2003). Each of the three species occurs at the focal mainland locations, as well as on Santa Cruz, Santa Rosa and Santa Catalina Islands. Only two of the three species, B. carinatus and N. pulchra, occur on Anacapa and San Miguel Islands (Junak et al. 1997). Anacapa Island is a chain of three islets; we sampled East Anacapa, which is most accessible. Habitats vary among these locations (Heady 1977). Native grasses are commonly located in patchy populations in coastal prairies and riparian areas of the Channel Islands and Vandenberg Air Force Base (VAFB), and in oak woodland savannahs at Sedgwick Reserve (Sedgwick).

Field collections

We located and georeferenced populations of each bunchgrass species during spring and summer of 2002 and 2003 with the goal of representing the geographic range of the three species across the islands and within the mainland study region (Fig. 1). Sampling was restricted to sites where one or more of the three species was present in order to sample species‐associated soil properties. At each site, three soil cores (7 cm × 20 cm) located between 10 and 20 m apart within the sampled population were collected using an AMS slide hammer (AMS, American Falls, Idaho, USA). Cores were bagged separately and transported to the University of California at Santa Barbara (UCSB) for subsequent analyses. We sampled a total of 118 sites, including twelve sites where two of three species overlapped (Table 1).

The total number of populations of each species sampled among locations. Altogether, there were 130 collections, of which two of three species overlapped at 12 sites. This resulted in 118 sites for analysis.

Soil sampling and analyses

Soil cores were sieved to 4 mm shortly after collection and stored at 4°C until the initiation of rewetting treatments approximately 60 days later. Water holding capacity (WHC) was determined by weighing a sample of 30 g of soil saturated with distilled water, and then drying the sample overnight at 75°C prior to measuring dry weight. Soils were subsequently adjusted to 50% WHC for nitrogen mineralization analysis, which was conducted for duplicate samples of 10 g of soil after the first day and also after a 14‐day incubation period at 20°C (Weintraub and Schimel 2003). We analyzed extractable ammonium (NH₄⁺) with the diffusion method, and nitrate (NO₃⁻) using the Griess‐Ilosvay reaction after reduction with cadmium (Lachat Instruments, Milwaukee, Wisconsin, USA; Lachat Methods 12‐107‐04‐1‐B and 31‐107‐06‐5‐A). Net nitrogen mineralization rates were calculated by subtracting the initial soil nitrogen (NH₄⁺ plus NO₃⁻) content from the sum of values for NH₄⁺ and NO₃⁻ measured in incubated soil.

Soil pH was measured for a 1:2 soil and water suspension that was shaken for 30 minutes and allowed to settle for five minutes. Samples were also analyzed for cation exchange capacity (CEC) including exchangeable calcium, magnesium, sodium and potassium (measured in milliequivalents per 100 g soil) using the NH₄‐OAc method buffered at pH 7.0 (Lavkulich 1981). Total carbon and nitrogen content (mg/L) was determined by analysis of dried samples on an elemental analyzer (Fissons Instruments, Milan, Italy) using an acetanilide standard. Duplicates were run every 10 samples to check accuracy of analysis (difference of 1–3%). Particle size analyses were performed at the Division of Agriculture and Natural Resources Analytical (DANR) Laboratory at the University of California Cooperative Extension in Davis, California following Sheldrick and Wang (1993).

We analyzed each of the three cores per site separately and subsequently averaged results for all variables with the exception of exchangeable cations and particle size. For those two variables, we combined equal portions of 100–200 g of each subsample to create a single sample for each site prior to analyses. Average soil pH was calculated as the mean of H+ concentration (pH = −log [H+]). Altogether, 15 soil variables were measured for use in subsequent data analyses. These variables included soil pH; WHC; total carbon and nitrogen; exchangeable calcium(X‐Ca), potassium(X‐K), magnesium(X‐Mg) and sodium (X‐Na); CEC; percent clay, sand and silt; initial NH₄⁺ and NO₃⁻; and net nitrogen mineralization (or delta N).

Data analysis

Soil variables were compared to species presence‐absence data using canonical redundancy analysis (RDA) and variation partitioning (Legendre and Legendre 1998). Prior to RDA, species values were Hellinger transformed to account for the large number of zeros in the dataset (Legendre and Gallagher 2001). To improve symmetry of the environmental data and meet expectations of statistical tests, the ratio of carbon and nitrogen, cations including the ratio of magnesium to calcium, and nitrogen mineralization values were ln(x + 1) transformed. Particle size data were arcsine square‐root transformed. Ordination methods were performed with R software version 2.15.1 (R Development Core Team 2009).

RDA was first conducted for the full dataset (island and mainland) using a matrix of spatial data created by recoding the five islands and two mainland locations as binary “dummy” variables to indicate the presence or absence of a sampled population at each island or mainland site (Legendre and Legendre 1998, Borcard et al. 2011). We chose this method due to the irregular sampling scheme required as a result of the variable geographic distances among the five islands and two mainland locations. Environmental data were altered to include standardized soil variables as well their second and third order polynomials to model nonlinear relationships in the dataset (Makarenkov and Legendre 2002). To reduce variance inflation, we selected among soil variables by forward selection with randomization tests (9,999 permutations) using the “packfor” R package (http://r‐forge.r‐project.org/R/?group_id=195) based on methods described by Blanchet et al. (2008). Exchangeable sodium was not included in RDA analysis of the full dataset due to multiple missing values for Sedgwick samples.

Once we determined the variables to include in the environmental (or soil) dataset, we partitioned variation among explanatory data and covariables using the varpart function in the “vegan” library (Oksanen et al. 2008). This method computes the proportion of variation in species abundance explained by: (1) the environmental dataset, (2) the spatial dataset, and (3) the fraction explained jointly by all spatial and environmental variables (Borcard et al. 1992, Peres‐Neto et al. 2006). As a final step, we ran the RDA using Hellinger‐transformed species response data, the reduced soil explanatory data, and spatial covariables (Borcard et al. 2011).

To better quantify the contribution of spatial and environmental variables to species abundance, we conducted RDA for a subset of samples that described the four islands in Channel Islands National Park (Anacapa, Santa Cruz, Santa Rosa and San Miguel Islands; Fig. 2). The geographical distance among the northern Channel Islands (CI) was small relative to other locations in the complete dataset, allowing study of finer spatial functions. Namely, the minimum spanning tree by which any two sites are connected using a single path had a truncation distance of 18.62 km. The truncation distance represents the single largest distance in the tree. Spatial coordinates among the CI sites were obtained with the geoXY() function in the SoDA package for R (Chambers 2008). Coordinates were used to construct a matrix of spatial explanatory variables using the method of Moran's eigenvector maps (MEM; Borcard and Legendre 2002, Dray et al. 2006). Distance‐based MEM is a method that identifies orthogonal spatial variables across multiple scales and can model both positive and negative correlations among species abundance data (Borcard et al. 2011). MEM spatial variables were only calculated for the four northern islands because these variables may be distorted if the truncation distance among sampled sites is too large (Borcard and Legendre 2002, Borcard et al. 2011).

Distribution of sampling points among the four Channel Islands for MEM analysis.

The CI‐only species data were Hellinger transformed, and standardized soil data included second and third order polynomials similar to analysis of the full dataset. Forward selection was conducted separately for spatial and soil variables. In the case of the soil dataset, we modified the analysis to use the significance value P ≤ 0.10 as the cutoff with the goal to identify variables that approached significance at P = 0.05. Variation partitioning was also conducted for the CI dataset. Final data for RDA analysis included the Hellinger‐transformed species data, the reduced soil environmental data, and matrices of (positive or negative) MEM spatial variables.

While forward selection reduced the number of explanatory variables in analyses, correlations among remaining variables may have inflated standard errors and obscured significant environmental and spatial patterns in the data. We therefore calculated variance inflation factors, which measure the degree of multicollinearity for each variable given the presence of other variables, for each RDA dataset. Inflation factors less than 10 indicated that correlations among variables were minimized (Borcard et al. 2011).

Our goals included testing the application of these methods, and so we examined the potential for individual soil variables identified by canonical ordination to accurately predict species presence or absence at a site. We conducted these analyses with logistic regression of selected soil variables separately for each of the three species. Logistic regression analysis of transformed soil variables and binary species data with a logit link function was performed in JMP version 9.0 (SAS Institute, Cary, North Carolina, USA).

Results

Island and mainland locations

Forward selection for samples among the five islands and two mainland locations included seven soil variables that explained a significant (0.05 level) portion of the variance in the distribution data. The first four variables were the percentage clay, magnesium‐calcium ratio, CEC and the net nitrogen mineralization rate. The remaining three selected variables were polynomial functions that included cubic values of percentage silt and the magnesium‐calcium ratio, and the quadratic value of CEC. Variance inflation factors were ≤2.66, and were well within the accepted range of values for multivariate analysis.

Selected soil variables explained 27% (adjusted R² or R_a²) of the variation in species presence‐absence data when modeled with spatial covariables. Spatial covariables, however, explained none of the variation in the species data. Permutation tests were significant for environmental variables (P = 0.0001) but indicated no effect of spatial variables (P = 0.9208). The first two canonical axes of the RDA were significant, and explained 19.0% of the variation in the bunchgrass distributions (Table 2, Fig. 3). Namely, the first axis described soil variables of CEC, the magnesium‐calcium ratio and % clay, and net nitrogen mineralization. The second canonical axis described polynomials of silt and CEC. The RDA biplot of results distinguished among the three native bunchgrass species (Fig. 3). Specifically, the presence of N. pulchra was correlated with clay‐rich, serpentine soils. Presence/absence of B. carinatus was associated with CEC, although the relationship was non‐linear. Elymus glaucus presence was associated with soil texture (silt), low magnesium‐calcium ratios, and net nitrogen mineralization.

RDA biplot for the three native bunchgrass species (B. carinatus: BC, E. glaucus: EG, and N. pulchra: NP) and seven selected soil variables. Soil variables are magnesium‐calcium ratio (Mg.Ca) and its cubic polynomial, cation exchange capacity (CEC) and its quadratic polynomial, percentage clay, and the cubic polynomial of percentage silt. Net nitrogen mineralization is represented by delN. Scaling = 1.

Results of analyses including permutation tests and canonical coefficients of the two significant axes identified by each RDA. Quadratic (second‐order) and cubic (third‐order) polynomials of soil variables are indicated.

Channel Islands National Park

Spatial data for analysis of the four northern Channel Islands included 70 MEM eigenfunctions. Forward selection retained 11 of the original 70 spatial variables, but none of these represented positive spatial autocorrelation (P = 0.272, forward selection of first positive MEM variable). The absence of significant positive spatial variables agreed with results for the complete island and mainland dataset, suggesting a lack of spatial structuring in these grasslands. Although the use of negative eigenfunctions in community analysis is not well‐defined, negative MEM variables may model biotic interactions such as competition and niche differentiation (Borcard et al. 2011). We therefore ran RDA for the forward selected soil variables and Hellinger‐transformed species data both with and without negative MEM spatial variables.

Forward selection of soil variables for samples among the four Channel Islands included the magnesium‐calcium ratio, percentage clay and CEC. Also included were quadratic values of initial NH₄⁺ and CEC. The value of alpha was 0.05 or less for each of these variables despite the 0.10 cutoff. RDA without spatial functions was highly significant (P = 0.0001) and soil variables explained 22% (R_a²) of the variation in species abundance data. Variance inflation factors were ≤1.59 for this analysis. The first two canonical axes were significant and explained 17.0% of the variation in the bunchgrass abundance data (Table 2). The first axis described CEC, the magnesium‐calcium ratio and % clay. The second canonical axis described polynomials of initial NH₄⁺ and CEC. RDA biplot results indicated separation between N. pulchra and B. carinatus, with weaker evidence for an association between soil characters and the presence/absence of E. glaucus (Fig. 4).

RDA biplot for the three native bunchgrass species and five selected soil variables for the four Channel Islands. Soil variables are the magnesium‐calcium ratio (Mg.Ca), cation exchange capacity (CEC) and its quadratic polynomial, initial NH4+, and percentage clay. Scaling = 1.

Each of the first two RDA analyses identified a significant association between soil properties and species presence‐absence data, but residual eigenvalues were high (>70%). The addition of negative MEM spatial functions to the model reduced residual eigenvalues. RDA was significant for soils (P = 0.0008), and variance inflation factors were low (≤1.91). Soil variables explained 9%, and spatial MEM functions 29%, of the variation in species abundance data. The variance fraction explained jointly by spatial and soil variables was 13%, and residual eigenvalues were reduced to 49%. In sum, the 22% of the variation in species data explained by soil properties in the CI analysis was divided into 9% explained solely by soil variables and 13% explained by spatially structured soil data. The first two canonical axes were significant, and N. pulchra populations were correlated with the percentage clay and magnesium‐calcium ratio (Table 2, Fig. 5). Bromus carinatus was correlated with initial NH₄⁺ and the quadratic value of CEC. There was some evidence for the association of E. glaucus and soil CEC (Fig. 6).

RDA fitted site scores of the first (A) and second (B) significant canonical axes in MEM analysis plotted as bubble maps. Square size is proportional to the absolute values (see legend) and positive values are identified by black squares while negative values are represented by white squares. RDA scores are a linear combination of explanatory variables included in analysis.

RDA biplot for the three native bunchgrass species and five selected soil variables for the analysis of four Channel Islands and 11 MEM spatial eigenfunctions. Soil variables are described in Fig. 4. Scaling = 1.

Logistic regression

Forward selected soil variables determined in RDA analysis were run in separate logistic regression analyses for each species and results indicated that individual soil variables can predict species presence or absence at a site. Results reported here list values derived from logistic likelihood ratio tests tables (chi‐square degrees of freedom = 1). The two variables with the greatest overall effect were the magnesium‐calcium ratio and nitrogen mineralization rate. There was also some effect of soil texture, which was specific to the percentage clay. Bromus carinatus populations were negatively correlated with magnesium‐calcium ratios (χ² = 6.87, P = 0.0088) and somewhat positively correlated with nitrogen mineralization rates (χ² = 3.72, P = 0.0539). Elymus glaucus distributions were negatively correlated with the percent clay (χ² = 3.90, P = 0.0484) and positively correlated with nitrogen mineralization rates (χ² = 5.79, P = 0.0161). The greatest number of variables predicted N. pulchra distributions (Mg/Ca, χ² = 4.99, P = 0.0254; N‐min, χ² = 9.19, P = 0.0024; %clay, χ² = 13.28, P = 0.0003). The presence of N. pulchra was positively correlated with magnesium‐calcium ratios and percent clay, and negatively correlated with nitrogen mineralization rates.

Discussion

Determining how to effectively reestablish California native grasses requires an understanding of the environmental factors that contribute to the persistence of focal species. However, given that large stands of native‐dominated grasslands are no longer extant, identifying those factors is difficult. Our goal for this study was to assess the relative contribution of soil properties and spatial structure to species distributions. Canonical analysis confirmed that soil properties explained a significant proportion of the variance in the distributions of the three species, but no positive spatial autocorrelation was detected. This indicated that contagious processes such as dispersal factors (e.g., wind and animals) did not explain species occurrence among sampled grasslands (Peters et al. 2006, Borcard et al. 2011). These three species are capable of dispersing to all sites among sampled island and mainland grassland locations, and environmental factors subsequently select for their persistence (Soininen et al. 2007, Kramer et al. 2011).

The absence of dispersal limitations among mainland and Channel Island sites appears to contradict a prior study of native grassland recruitment, which concluded that low rates of bunchgrass seed production and restricted dispersal affected the competitiveness of these species in invaded grasslands (Seabloom et al. 2003). The relative importance of spatial and environmental factors, however, is dependent on scale and our work focused on existing native grasslands (Legendre and Legendre 1998). We evaluated soil characteristics associated with species presence and absence at sites currently occupied by remnant native bunchgass populations. Consequently, these results do not address the extent to which the three species can disperse to unoccupied (non‐native) sites, particularly if species are niche‐limited (Moore and Elmendorf 2006). Moreover, recruitment limitation noted in earlier studies of California grasslands may be caused by low rates of seed production, rather than poor rates of dispersal among existing seed pools (Moore et al. 2011).

The presence of significant negative spatial autocorrelation among sites for each species may be a consequence of species ecology, or a sampling effect. We sampled grassland sites specific to these three species in order to characterize environmental variables relevant for their persistence. As a result, many sites represented one species out of three, with 12 of 118 collections representing more than one native bunchgrass population. If sampling did not identify all species at a site, the sampling protocol may represent a source of experimental error. However, sites were revisited frequently and occurrence data were corroborated during each visit. It is therefore likely that data reflect tangible spatial and biotic processes, perhaps due to niche differentiation among species (Leibold and McPeek 2006, Aanderud and Bledsoe 2009). This explanation is supported by the significant environmental variables identified in subsequent ordination and logistic regression analyses. For example, the distributions of N. pulchra and B. carinatus were negatively correlated among sites for soil magnesium‐calcium ratio and nitrogen mineralization potential. In effect, the two species occupied contrasting environments.

Although we identified edaphic factors that predict species presence or absence, this does not confirm causation. Plants alter their environment by redistributing soil nutrients via root structures and litter accumulation as well as by producing allelopathic compounds (Weidenhamer and Callaway 2010, Lankau et al. 2011). In the case of abiotic soil properties, mineral and physical factors are a recognized cause of non‐random species distributions (Brady and Weil 2002). The question remains, however, whether biological and chemical soil properties such as nitrogen mineralization and cation exchange capacity represent causes or consequences of plant species distributions. Sites occupied by B. carinatus and E. glaucus were often more mesic and situated on deeper soils relative to sites occupied by N. pulchra (K. Hufford, personal observation). In addition, B. carinatus and E. glaucus did not grow as monocultures and instead were patchily distributed. These characteristics predict that biological and chemical soil properties are a function of the larger community, and may also be a causal factor of bunchgrass distribution in a nutrient‐poor environment (John et al. 2007).

In the case of the four northern Channel Islands, 13% of the 22% of the variance explained by soil variables was spatially distributed. Spatial heterogeneity may result from patchy distributions of soil properties, such as the occurrence of serpentine outcrops among sedimentary soils (Borcard et al. 2011). Serpentine soil is widely known to drive species distributions (Kruckeberg 1951, Schechter and Bruns 2008) and the ratio of magnesium to calcium was a significant factor determining species presence or absence in both ordination and logistic regression analyses. Bromus carinatus was less likely to occur on serpentine soils while N. pulchra populations appeared to favor serpentine substrate. Nassella pulchra was once considered to be a historically dominant species in California grasslands (Heady 1977), but later research noted that N. pulchra is adapted to disturbance, suggesting it is more likely a marginal species (Bartolome and Gemmill 1981, Lombardo et al. 2007). Invasive species may have forced N. pulchra populations to occupy serpentine soils as a refuge from intense competition for resources (Harrison et al. 2003). Results were not as conclusive for E. glaucus, possibly reflecting lower sample sizes and the need for additional environmental data as its distribution may be influenced by biotic or abiotic factors not included in this study. We note that island populations of E. glaucus were collected after extensive searches and were representative of the observed distribution.

The soil properties identified in canonical ordination closely matched logistic regression results, but were not identical in all analyses. This outcome may indicate the uneven predictive role of soil properties across a species range, as well as the varying power of different statistics. In addition, while environmental data explained a significant proportion (22–27%) of the variance in species abundance, a majority of the variation in the dataset was unexplained. Our set of environmental variables likely missed measures such as phosphorus and moisture content, as well as biotic interactions such as competition and plant‐microbial relationships, that contribute to species distributions (Seabloom et al. 2003, Hawkes et al. 2005, Hausmann and Hawkes 2009). As a result, soil properties play a role in species distributions in this landscape, but other factors will likely also affect establishment. A next step would be to verify the ability of edaphic factors to predict revegetation success, supporting the use of individual soil variables in restoration planning.

In sum, the significant associations between our focal species' distributions and magnesium‐calcium ratios, percentage clay and nitrogen mineralization potential are consistent among analyses, and suggest that differences among the preferred soil habitats for the three species can assist in restoration of native bunchgrasses. This result offers potential tools for restoration practitioners who could incorporate knowledge of the association of key species with soil properties into decision making for seed mixes and seeding plans. A combination of additional research and practical testing would be useful to determine if the variance in species' occurrence explained by soil properties is sufficient to influence restoration outcomes. Ultimately, practitioners will benefit from an understanding of local and regional patterns of diversity in contemporary, degraded ecosystems, and this understanding may well assist species selection for each restoration site.

Acknowledgments

We thank Allen Doyle and Tony Hartshorn for assistance with soil analyses, and Oliver Chadwick for the generous use of laboratory facilities. Field collections were completed with the assistance of Sarah Chaney and Dirk Rodrigues at Channel Islands National Park, Randy Miller at VAFB, and Mike Williams at UC Sedgwick Reserve. The Catalina Island Conservancy provided extensive field support. Numerous undergraduates assisted with sample processing, including Laura Dell'Olio, Drew Madrigal, and Katherine Davis. Analyses were assisted by a Master Class led by Pierre Legendre in Perth, Australia. This research was supported by the National Parks Ecological Research Fellowship Program, a program funded by the National Park Foundation through a generous grant from the Andrew W. Mellon foundation.