Non-Random Assembly of Bacterioplankton Communities in the Subtropical North Pacific Ocean

The exploration of bacterial diversity in the global ocean has revealed new taxa and previously unrecognized metabolic potential; however, our understanding of what regulates this diversity is limited. Using terminal restriction fragment length polymorphism (T-RFLP) data from bacterial small-subunit ribosomal RNA genes we show that, independent of depth and time, a large fraction of bacterioplankton co-occurrence patterns are non-random in the oligotrophic North Pacific subtropical gyre (NPSG). Pair-wise correlations of all identified operational taxonomic units (OTUs) revealed a high degree of significance, with 6.6% of the pair-wise co-occurrences being negatively correlated and 20.7% of them being positive. The most abundant OTUs, putatively identified as Prochlorococcus, SAR11, and SAR116 bacteria, were among the most correlated OTUs. As expected, bacterial community composition lacked statistically significant patterns of seasonality in the mostly stratified water column except in a few depth horizons of the sunlit surface waters, with higher frequency variations in community structure apparently related to populations associated with the deep chlorophyll maximum. Communities were structured vertically into epipelagic, mesopelagic, and bathypelagic populations. Permutation-based statistical analyses of T-RFLP data and their corresponding metadata revealed a broad range of putative environmental drivers controlling bacterioplankton community composition in the NPSG, including concentrations of inorganic nutrients and phytoplankton pigments. Together, our results suggest that deterministic forces such as environmental filtering and interactions among taxa determine bacterioplankton community patterns, and consequently affect ecosystem functions in the NPSG.


INTRODUCTION
The bulk of the world's ocean is oligotrophic and, in terms of both biomass and ecosystem function, is dominated by microorganisms (Cho and Azam, 1988). Microbes with diverse, yet specific physiological requirements and ecological roles mediate nearly all energy and material fluxes within this vast biome (Azam, 1998;Herndl et al., 2008). These fluxes have been suggested to be products of multiple interactions within the microbial community, rather than those of single populations (Little et al., 2008;Strom, 2008). Therefore, an understanding of the processes that determine the composition of bacterioplankton communities is imperative in order to formulate accurate estimates and models of open ocean biogeochemistry. Despite increasing efforts to understand emergent patterns in the abundances and distributions of microbial taxa in seawater (e.g., Morris et al., 2002;Fuhrman et al., 2006;Johnson et al., 2006;Carlson et al., 2009;Eiler et al., 2009), the mechanisms that determine patterns observed in entire microbial communities remain largely unknown (Fuhrman and Steele, 2008;Steele et al., 2011). This stands in contrast to macroorganismal ecology, where the ecological mechanisms that underlie and maintain community patterns (beta-diversity) and (alpha-) diversity (often termed "assembly rules") have been a central theme (e.g., Keddy and Weiher, 2001). Patterns in alpha-and beta-diversity of macroorganisms have been suggested to be the result of the interplay between regional and local processes: regional processes include dispersal, historical legacies, speciation, extinction, and geographic and climatic conditions, whereas local processes include interactions among organisms, and between organisms and their local environment.
Recently, microbial ecologists have begun to adopt many theoretical ecological frameworks upon which the microbial world can be explored and, eventually, predicted (e.g., Horner-Devine et al., 2004;Green and Bohannan, 2006;Martiny et al., 2006;Prosser et al., 2007). For example, microbial ecologists are now exploring whether microbial communities are distinct in (or specific to) different habitats, and if metrics of microbial diversity (alpha-and beta-diversity) show explainable patterns comparable to those of macroorganisms (e.g., Fuhrman et al., 2006;Horner-Devine et al., 2007;Prosser et al., 2007).
The balance between stochastic and deterministic processes in determining the structure of pelagic microbial communities in the open ocean has received little direct attention. If stochastic processes dominate, then open ocean communities may appear as random combinations of taxa, whereas, if deterministic processes dominate, then predictable patterns of taxa distributions and abundances may be observed. In order to test if community assembly is mainly established non-randomly and thus by deterministic forces, computer simulations can be used that compare simulated, randomly assembled communities with experimental data. Such an assessment can be achieved by analyzing taxa co-occurrence patterns via the approach proposed by Diamond (1975) and modified by Gotelli and McCabe (2002). Using this approach, non-random co-occurrence patterns have been shown to be the rule for macroorganisms (Gotelli and McCabe, 2002), and have also been demonstrated for microorganisms in a number of habitats (Horner-Devine et al., 2007). Environmental filtering and competition (Diamond, 1975;Gotelli and McCabe, 2002), and/or biogeographic and evolutionary history (Vuilleumier and Simberloff, 1980;Cracraft, 1988) are thought to be important deterministic processes that determine community composition.
However, in order to explore the mechanisms that determine community composition, it is not sufficient to describe patterns in community structure and to demonstrate that co-occurrence patterns are random or non-random; the explicit mechanisms that underlie the emergent community patterns need to be identified. Statistical tools allow for the exploration of community drivers through the correlation and/or correspondence among taxa and between taxa and environmental variables.
In the present study, we aimed to explore patterns in bacterioplankton community structure within water samples collected from various depths over a 3-year period at Station ALOHA in the North Pacific subtropical gyre (NPSG), the world's largest oligotrophic marine ecosystem (Karl, 1999). Our analysis revealed that deterministic mechanisms prevail in assembling the bacterioplankton community at Station ALOHA, and that associations among taxa, as well as between taxa and environmental variables such as nutrient and pigment concentrations, are common and predictable.

SITE DESCRIPTION AND ANCILLARY DATA
Station ALOHA (22˚45'N, 158˚00'W) is located in the NPSG 100 km north of Oahu, Hawaii. The Hawaii Ocean Time-series (HOT) program 1 (Karl and Lukas, 1996) has made repeated observations of the hydrography, chemistry, and biology of the water column at Station ALOHA approximately once per month since October 1988. Sampling and analytical protocols are described elsewhere Karl and Lukas, 1996;Letelier et al., 1996;Karl et al., 2001;Church et al., 2005). Metadata used for multivariate statistical analyses were downloaded from the HOT Data Organization and Graphical System (HOT-DOGS) website 2 . Environmental variables used in this study are listed in Table A2 in Appendix.

TERMINAL RESTRICTION FRAGMENT LENGTH POLYMORPHISM ANALYSIS
Two microliter of each DNA extract (equivalent to ca. 1-100 ng of genomic DNA) was used as template for the polymerase chain reaction (PCR) amplification of bacterial small subunit (SSU) ribosomal RNA (rRNA) genes. Individual 50 μL reactions contained 1× PicoMaxx reaction buffer (Stratagene, La Jolla, CA, USA), 200 μM of each deoxynucleotide, 200 nM each of the fluorescently labeled general bacterial SSU rRNA gene oligonucleotide primer 27F-B-6FAM (5 -[6∼FAM]AGRGTTYGATYMTGGCTCAG-3 ) and universal SSU rRNA gene oligonucleotide primer 1492R (5 -GGYTACCTTGTTACGACTT-3 ; Lane, 1991) and 0.625 units PicoMaxx thermostable DNA polymerase (Stratagene). Cycling conditions consisted of an initial denaturation step of 95˚C for 5 min, followed by 35 cycles of 95˚C denaturation for 30 sec, 55˚C annealing for 1 min, and 72˚C extension for 2 min, with a final extension step at 72˚C for 20 min. PCR products were purified using the QIAquick Multiwell PCR Purification system (Qiagen Inc.), and digested in 20 μL reactions containing 4 μg acetylated bovine serum albumin (BSA), 1X MULTI-CORE enzymatic reaction buffer (Promega, Madison, WI, USA), and 10 units of HaeIII restriction endonuclease (Promega). After a 6-h incubation at 37˚C, the digests were terminated by incubation at 70˚C for 20 min, and purified using the Millipore MultiScreen Assay System (Millipore Corp.) in conjunction with Sephadex G-50 Superfine (GE Healthcare, Piscataway, NJ, USA). Approximately 100 ng of the purified products were electrophoresed on an ABI 3100 Genetic Analyzer (Applied Biosystems, Foster City, CA, USA). Operational Frontiers in Microbiology | Aquatic Microbiology taxonomic units (OTUs) were identified as peaks between 30 and 600 bp in length, using GeneMapper v4.0 software (Applied Biosystems, Foster City, CA, USA). Terminal restriction fragment length polymorphism (T-RFLP) profiles were first run through an in-house MATLAB (The MathWorks, Natick, MA, USA) script that corrects for instrumental drift while rounding fractional peaks to whole numbers. In brief, this script identified peaks with ambiguous decimal increments (i.e., x.3-x.6) and ensured that they were rounded to the same integer throughout all profiles by using a majority rule (i.e., rounded up if most were > x.5 and vice versa). Non-ambiguous decimal increments (i.e., <x.3 and >x.6) were simply rounded to their nearest integer. The profiles then underwent a limited manual curation in order to further compensate for instrumental drift, and finally normalized using a variable threshold (Osborne et al., 2006). Peaks were transformed into relative units by dividing integrated peak areas by the total peak area for that sample. Diversity (evenness) was estimated from the relative units using the Gini coefficient (Wittebolle et al., 2009). The Gini coefficient provides a proxy for the equality of taxa distributions: the higher the Gini coefficient, the more uneven a community is.
Putative identification of abundant OTUs was performed using clone library data obtained from both coastal and open ocean sites in the vicinity of Hawaii, including Station ALOHA. Representative clones were characterized by Sanger sequencing and fragment sizing after restriction digests to correlate phylogenetic identity with terminal restriction fragment length. Sequences were deposited in Genebank under accession number JN166096 to JN166367.
Statistical analyses required the exclusion of samples with missing values. In order to account for such gaps in the availability of metadata, environmental variables were grouped into biological, chemical, and physical data classes in order to increase sample size for analyses ( Table A2 in Appendix). Within PRIMER 6, BIOENV procedures were used to explore the correspondence between microbial community structure and environmental metadata (Clarke and Warwick, 2001). In this analysis, the environmental variables were first Log(X + 1) transformed and normalized by subtracting the mean from each entry of a single variable and dividing by the standard deviation for that variable. This allowed the identification of variables or combinations of variables that best corresponded with the community data. As this analysis was based on a limited number of samples due to missing data, the best fitting combinations were re-analyzed after adding samples with non-missing values (in the identified variable combination) using the RELATE test (Clarke and Warwick, 2001) in PRIMER 6 (PRIMER-E Ltd.). The Analysis of Similarity (ANOSIM in PRIMER 6; Clarke, 1993;Clarke and Warwick, 2001) test was used to examine for differences between bacterial communities from different depth horizons in the water column.
Two taxa co-occurrence indices were calculated using the EcoSim 700 software package (Gotelli and McCabe, 2002). The Checkboard index (C-board) calculates the number of taxa that never co-occur, which can be interpreted as evidence for competitive exclusion (segregation) or differences in habitat preference. The second index (C-score) is calculated by: where Ri and Rj are the number of sites where taxon i and j occur, and S is the number of sites containing both taxa. The C-score is then calculated from the mean number of all pair-wise CUs. Monte Carlo simulations were used to test if the observed co-occurrence measures were random following the most conservative permutation procedure. Following Gotelli and McCabe (2002), we calculated a standardized effect size (SES) that allowed comparing the degree of co-occurrence across temporal and spatial scales.
Spearman's rank correlation coefficient (ρs) was calculated in MATLAB in a pair-wise fashion between OTUs, using the relative contribution of each OTU to the total community for each sample. Only OTUs that were detected in more than two independent samples were included in the analysis (n = 397). A significance level of p < 0.05 was used to calculate the percentage of significant positive and negative correlations from the total number of pairwise comparisons. In addition, Bonferroni correction was used to address the problem of multiple comparisons. In the same way, ρs between relative contributions of each OTU and depth were calculated.
To test if the number and distribution of significant pair-wise correlations was different from random, Monte Carlo simulations were performed. Following a permutation procedure where the rank of OTUs within samples was randomized, a total of 1000 randomized matrices were created and subsequently used in pair-wise correlations. The maximum and minimum ρs values of correlations within each 0.01 interval (from 1 to −1) were determined from all 1000 randomized matrices.

CO-OCCURRENCE PATTERNS AND CORRELATIONS BETWEEN OTUs
A total of 461 different OTUs, each represented by a specific terminal restriction fragment length, were detected in 296 samples collected over 3.5 years (November 2004 to February 2008) from depths ranging from 5 to 4000 m at Station ALOHA. Of these, 397 OTUs were present in at least two T-RFLP profiles, and so were used for further analyses. Co-occurrence indices calculated from this dataset differed significantly from random simulations: more OTU pairs exhibited larger C-scores and C-boards than expected by chance ( Table 1). A higher C-score than expected by chance indicates low co-occurrence, whereas a lower C-score than expected by chance indicates high cooccurrence. Low co-occurrence can be interpreted as evidence for competitive exclusion (segregation), or differences in habitat preference. Of the Spearman's rank correlation coefficients (ρs) from all pair-wise comparisons of OTUs (n = 78,607), 27.3% were significant at a 5% significance level, of which 6.6% were negatively correlated (ρs < 0) and 20.7% were positively correlated (ρs > 0; Figure 1A). Both significantly positive and negative correlations occurred more often than expected by chance when compared to simulations where the rank order of species were randomized within a community. In addition, maximum ρs values were much higher for the observed matrix compared to the simulated matrices ( Figure 1A). Using Bonferroni correction for addressing the problem of multiple comparisons, the number of significant correlations decreased from 21,485 to 3,276 (3.5% significant positive and 0.5% significant negative correlations).

Frontiers in Microbiology | Aquatic Microbiology
On average, 108 significant correlations were observed per OTU; 22 OTUs were correlated with more than 50% of the 397 OTUs included in the dataset. This subgroup of 22 included the most abundant OTUs, which are putatively identified as Prochlorococcus and SAR11 subclades IB and IA. Other abundant OTUs putatively identified as SAR116 subclade II and SAR11 subclade IIB were significantly correlated with >150 OTUs. However, many OTUs contributing <1% in relative abundance (on average) to the total bacterioplankton community exhibited a high number of significantly correlations with other OTUs in our dataset ( Figure 1B). Hence, we did not observe a significant correlation between the average relative abundance of an OTU and the number of correlations it exhibited.
Variation in the C-score along the depth gradient was investigated in detail by calculating the SES from C-score calculations grouped by depth (SEScscore; Figure 2). A positive SEScscore significantly different from random (p < 0.01) was observed at all depths, indicating low co-occurrence. The SEScscore was not correlated with depth, but exhibited intriguing variability along the vertical axis (Figure 2). The highest value was estimated at 400 m, and decreased again toward the oxygen minimum zone at ca. 800 m. The SEScscore did not vary with the size of the data set, as measured by either the number of samples or the number of different OTUs included in a given incidence matrix (ρ = 0.103, p = 0.696, and ρ = 0.359, p = 0.137, respectively). It could also be argued that the number of OTUs represents an estimate of FIGURE 2 | (A) Depth-specific distribution of oxygen and chlorophyll a concentrations from HOT cruises 165-200 (November 2004to February 2008. (B) Standardized effective size of C-score (SEScscore), average Bray-Curtis similarity from all pariwise comparisons within a single depth, and Gini coefficients from depth-specific communities. Average Bray-Curtis similarity and Gini coefficients were calculated from T-RFLP-based community structure profiles pooled by depth. SESscore was calculated by adapting the procedure proposed by Gotelli and McCabe (2002), using all time points for a specific depth.
www.frontiersin.org community richness. Hence, this apparent lack of relationship between richness and SEScscore, together with evenness and SEScscore also not being related (ρ = −0.279, p = 0.276) indicates that co-occurrence patterns are not related to alpha diversity at least in our study.
The Gini coefficient, an estimate of community evenness, ranged from 0.007 to 0.067 in the 296-sample dataset. These low values indicate a high degree of equality, as a value of 0 indicates complete equality while a value of 1 indicates maximal inequality among OTUs. Spearman's rank correlation analysis revealed a significant negative relationship between the Gini coefficient and depth (ρ = −0.27, p < 0.001) indicating an increase in community evenness with depth (Figure 2).

CORRESPONDENCE BETWEEN MICROBIAL COMMUNITY STRUCTURE AND ENVIRONMENTAL VARIABLES
In order to visualize potential spatial and temporal patterns in community structure, communities were plotted using NMS to two dimensions based on a Bray-Curtis distance matrix (Figure 3). ANOSIM and RELATE analyses corroborated initial visual observations of depth-specific clustering of bacterioplankton communities at Station ALOHA (R = 0.613, p < 0.001, and ρ = 0.67, p < 0.001, respectively). An R-value of 0 resulting from the ANOSIM analysis indicates that there were no differences among the communities from different depths, while an R-value of close to 1 would indicate that the dissimilarities between communities from different depths were larger than those from each individual depth. However, pair-wise ANOSIM revealed that community structure from different depth horizons in the upper 100 m were not significantly different whereas, below 100 m, all depth horizons were significantly different from each other ( Table A1 in Appendix). Exceptions were adjunct depths where, in most cases, differences were not significant ( Table A1 in Appendix). Average Bray-Curtis dissimilarities between communities from different depths increased with increasing distance along the vertical (i.e., depth) profile. The significant ρ value of the RELATE analysis indicated that communities from more distant depth horizons were more different than communities from adjunct depths. When Bray-Curtis similarities were averaged over the sampling period for each depth, values stayed rather constant above the mixed layer (Figure 2). In the depth layers where the mixed layer and chlorophyll maximum temporally oscillate, Bray-Curtis similarity averages were low. The relative abundance of nearly early half (46%) of all OTUs was significantly correlated with depth (p < 0.05), of which 32% were positive and 14% negative (Figure 4).
Unlike the clear depth-specific zonation of bacterioplankton communities, temporal patterns were not readily apparent in NMS plots (data not shown). However, we were able to statistically test for temporal patterns in bacterioplankton community structure, since the sampling schema spanned over more than three consecutive years. When months and seasons were used as discriminating factors in ANOSIM analyses, and when months, seasons, and Julian days were used in cyclic RELATE analyses, significant temporal patterns were observed at various depths. However, significant annual cycles in bacterioplankton community composition were only found in the upper 200 m of the water column ( Table 2) and were not detected in the deeper horizons.
BIOENV procedures were employed to search for relationships (correspondence) between community and environmental variability, and to identify the strongest relationships. Since the environmental dataset was incomplete, variables were divided into three categories: physiochemical, biogeochemical, and biological ( Table A2 in Appendix). Data from each category was used separately in the BIOENV procedure, and the variable combinations that best matched the patterns in community composition were identified and then analyzed by RELATE after inflating the number of samples (Table 3). Parameters from all three general categories revealed significant correspondence along the sampled spatial and temporal scales ( Table 3). This included depth, which represented the single variable that corresponded the best with community composition. Variation in bacterioplankton community structure in the upper 200 m of the water column was highly correlated with concentrations in heterotrophic prokaryotes, Prochlorococcus, Synechococcus, and pigments such as chlorophyll a and phaeopigments. In addition, inorganic nutrients such as nitrite, nitrate, and particulate carbon corresponded well with community composition. Most of these variables were also correlated with depth ( Table A3 in Appendix), and annual periodicity in these variables was either weak or not present (data not shown).

DISCUSSION
In this study, we sought to understand the mechanisms that underlie the patterns of community composition and maintain biodiversity over space and time (in other words,"to formulate community assembly rules"). Our results suggest that non-random patterns of co-occurrence are common between bacterioplankton lineages in the NPSG. This is consistent with observations made previously in other ocean provinces Gilbert et al.,   2009; Steele et al., 2011), and demonstrates that both the emergent depth-dependent and temporal patterns in bacterioplankton taxa distributions in the NPSG derive mainly from deterministic forces, like controlling environmental conditions and taxa interdependencies. Documenting non-random patterns of co-occurrence for members of natural microbial communities is an important www.frontiersin.org step in establishing generalities for community assembly, and confirms previous predictions regarding the non-random assembly of microbial communities in other environments (Horner-Devine et al., 2007). In this meta-analysis, no clear trends in co-occurrence patterns could be observed in relation to microbial habitats or phylogenetic groups but trends could be attributed to the variability in sample sizes and analysis methods of the different studies. In our analysis, SEScscore values were not related to sample size, even though the number of samples analyzed from each depth ranged from 6 to 33. In addition, proxies for alpha-diversity (GINI index and the total number of OTUs) and beta-diversity (Bray-Curtis similarity) were not related to SEScscore values, suggesting that community properties are not likely effecting cooccurrence estimates. Still, we observed depth-specific variability in co-occurrence indices (SEScscore) within the predominantly stratified water column of the NPSG that could be attributed to the presence of complex and depth-specific ecological forces structuring bacterioplankton communities in the oligotrophic ocean. High SEScscore values appear to coincide with steep, depth-specific gradients in such environmental variables as chloropigments and oxygen concentrations (Figure 2). When coupled with the observation that bacterioplankton communities are stratified in the NPSG, it appears that deterministic mechanisms, such as environmental selection and competition, are stronger in these deep horizons of the water column than, for example, stochastic mixing events that transport nutrients and organisms through the water column. In addition, these stratified bacterial communities were maintained consistently for 3.5 years. An exception to this was observed above the depth of the mixed layer, where communities were fairly homogeneous (Figure 3; Table A2 in Appendix). This stands in contrast to the deep mixing events that disrupt the stratified bacterioplankton communities at the Bermuda Atlantic Time-series Study (BATS) site in winter, producing large shifts in community composition in the epipelagic and upper mesopelagic (Morris et al., 2005;Treusch et al., 2009).

Frontiers in Microbiology | Aquatic Microbiology
While simple randomization tests like those applied in this study are useful statistical tools for recognizing non-random patterns of distribution (Gotelli, 2001), they do not allow for the identification of the causal underlying mechanisms shaping bacterioplankton communities in this system (Roughgarden, 1983). If the goal is to predict bacterioplankton community composition in the oligotrophic ocean, the causal underlying mechanisms need to be disentangled and their relative importance quantified. The application of statistical tools to identify relationships between community structure and abiotic and biotic environmental properties, as well as relationships among bacterial groups, is a first step toward such predictions (e.g., Fuhrman et al., 2006;Gilbert et al., 2009;Treusch et al., 2009). Using permutation-based procedures such as ANOSIM, we observed highly structured bacterioplankton communities along depth gradients, and a high number of individual OTUs either positively or negatively correlated with depth, allowing predictions on their distribution along the stratified water column. Numerous studies have demonstrated that specific groups of planktonic microorganisms are vertically stratified in the open ocean (e.g., DeLong, 1992;Giovannoni et al., 1996;Field et al., 1997;Church et al., 2005;DeLong et al., 2006;Johnson et al., 2006;Pham et al., 2008). While not surprising when considering the existence of strong environmental gradients along the vertical profile (Table A3 in Appendix), disentangling depth from correlations with key environmental parameters remains an important but unresolved task.
In addition to identifying that different depths harbor distinct bacterioplankton communities, we were able to characterize these differences in evenness, SEScscore and variance in community similarity ( Figure 2B). Community evenness increased with depth, whereas the variability in community composition (as indicated by Bray-Curtis similarities) decreased with depth, clearly showing that the epipelagic bacterioplankton community was more dynamic than the bathylopelagic community. Higher variability in epipelagic communities was also reflected by the subtle but statistically significant seasonal oscillations, suggesting re-occurring temporal patterns as proposed for other ocean sites Treusch et al., 2009). There was also a clear correspondence between epipelagic bacterioplankton community structure and environmental variables that do exhibit some seasonality at Station ALOHA (Karl, 2002). The relationship between bacterioplankton community composition and nitrite, nitrate, and particulate carbon concentrations and the abundance of Prochlorococcus, Synechococcus, and pigments indicate that resource availability (or bottom up control) plays a dominant role in the assembly of bacterioplankton communities in the NPSG (Torsvik et al., 2002). Previous studies also Frontiers in Microbiology | Aquatic Microbiology point to the importance of top down controls such as viruses and grazers (Thingstad and Lignell, 1997;Fuhrman, 1999;Hewson et al., 2006). Besides trophic interactions such as bacterivory and viral lysis, other forms of interdependencies have attracted little attention, especially antagonism and mutualism among open ocean bacterioplankton taxa.
Correlations among taxa may provide preliminary indications regarding the nature of bacterioplankton interactions. For example, a positive correlation may reflect syntrophy or other interdependency between two taxa, though it may also reflect similar responses to controlling environmental conditions, leading to the two taxa following a similar environmental trend. A negative correlation, on the other hand, may reflect either direct competition, or differences in growth constraints between taxa. Still, the higher number of positive than negative correlations among taxa observed in this study could partially result from a higher degree of mutualistic than antagonistic interactions. Metagenomic sequence data obtained along a depth profile in the NPSG support this view, as many genes with potential interactive functions have been identified, such genes encoding for chemotaxis, antibiotics, pili, and flagellate synthesis (DeLong et al., 2006). However, metabolic interactions between microorganisms resulting from the biotransformation of metabolites (e.g., from reduced to oxidized states and vice versa) do not rely on close spatial association, and hence may be far more important than interactions requiring close physical proximity. Many transformations in the carbon, nitrogen, and sulfur cycles are reliant on appropriate cellular and ambient microenvironmental redox conditions (see for example, Paerl and Pinckney, 1996, and references therein). Genome analyses have also emphasized the importance of metabolites for interactions in pelagic marine environments (Tripp et al., 2008). Although many of the identified correlations among taxon-pairs may be the result of such metabolic interactions, experimental validation is required in order to explore the causality behind the observed correlations. A logical starting point for such studies are the highly correlating taxa, which may represent groups with extensive and strong interactions to other taxa. The high relative abundance of some of these taxa (i.e., Prochlorococcus and SAR11) already suggests a key role in open ocean bacterioplankton communities.

CONCLUSION
The non-random temporal and spatial patterns in microbial community structure we observed are an important starting point for investigating microbial ecosystem services in the NPSG . For example, if interdependencies including mutualism, antagonism, trophodynamics (predation and parasitism), and competition determine community structure, then it can be speculated that these interactions directly underlay ecosystem processes by modifying pathways of energy and material flow, or indirectly by modifying the abundances of organisms or traits with strong ecosystem effects. Thus, understanding community interactions is essential if we want to predict ecosystem functioning of the oligotrophic ocean in an era of anthropogenic induced environmental change. www.frontiersin.org