Geometrical frustration of phase-separated domains in Coscinodiscus diatom frustules.

Diatoms are single-celled organisms with a cell wall made of silica, called the frustule. Even though their elaborate patterns have fascinated scientists for years, little is known about the biological and physical mechanisms underlying their organization. In this work, we take a top-down approach and examine the micrometer-scale organization of diatoms from the Coscinodiscus family. We find two competing tendencies of organization, which appear to be controlled by distinct biological pathways. On one hand, micrometer-scale pores organize locally on a triangular lattice. On the other hand, lattice vectors tend to point globally toward a center of symmetry. This competition results in a frustrated triangular lattice, populated with geometrically necessary defects whose density increases near the center.

Phase separation is widely appreciated as an important mechanism used by cells to control biochemical pathways (1, 2). However, cells also exploit phase separation to form functional structures with remarkable order and symmetry. For example, such structures can be found in the keratinaceous color-producing nanostructures of bird feathers (3, 4), the topography of pollen’s cellulosic cell wall (5), and the porosity of silicaceous diatom frustules (6–8). Some remarkable ordered porous structures can be found at a larger scale in sea urchin spines (9) and knobby starfish (10).

Diatoms are single-celled organisms with a beautifully structured cell wall made of silica, called the frustule. Their elaborate patterns have fascinated scientists and artists for hundreds of years (11). The diatom builds its wall by precipitating silica in a specialized organelle called the silica deposition vesicle (SDV) (12), templated by the phase separation of biopolymers, such as long-chain polyamines (LCPAs) (7, 13). The resulting structure is multiscale (7) and displays a remarkably high specific strength (14). The intricately patterned frustule is thought to offer protection from predators (15) and to control the transport of mass (16) and possibly, light (17). While much progress has been made in understanding the chemical processes underlying silica deposition, little is known about the frustule’s organization at larger length scales. Experiments in vitro have shown that appropriate mixtures of LCPAs, silaffins, and other proteins can spontaneously form a porous silica structure (18, 19). However, these structures are highly disorganized and show none of the long-range organization that characterizes frustules. The organization of these structures has been hypothesized to emerge spontaneously either through the self-assembly of phase-separated LCPA-rich domains (7) or through templating by the cytoskeleton (20, 21).

Here, we quantify the structural arrangement of micrometer-scale pores in valves of the diatom Coscinodiscus granii. We find that the pores locally form a triangular lattice. Globally, however, the structure is radially aligned. Other diatoms from the Coscinodiscus family also tend to form a local triangular lattice but with varying degrees of radial alignment. Valves of C. granii grown in the presence of germanium maintain radial alignment in the absence of significant hexagonal order. Local hexagonal order is consistent with the self-assembly of mobile uniformly sized domains with steric repulsions. The mechanism underlying radial alignment is unknown, but it is apparently modulated across species. Using a minimal Brownian dynamics simulation with hard disks and bond rotation, we created triangular lattices with various degrees of radial alignment. The simulated patterns span the range of morphologies observed in the Coscinodiscus family. While the competing tendencies of alignment and crystallization can readily be accommodated at the periphery of the frustule, geometric frustration leads to a highly disorganized core near the center.

Results and Discussion

We chose C. granii because its relatively large size and microscale pores enable optical microscopy for quantitative image analysis. In addition, it is one of the more established model diatoms (22–24). To characterize the structure of valves, the flat regions on either end of the frustule, we imaged their micrometer-scale features using spinning-disk confocal microscopy (100× numerical aperture 1.45). Starting with three-dimensional (3D) image stacks of the extracted rhodamine-labeled diatom valves (Fig. 1), we calculated the 3D centroids of the micrometer-scale pores in the valve (Fig. 1). From the side view, it is clear that the valve is mostly flat, with some curvature near the edge. As locating pores was less reliable in these areas (shown in gray), we excluded them from further analysis.

Fig. 1.

C. granii shows an evenly spaced pattern of micropores. (A) An average intensity projection of a confocal image stack for a representative valve. (Scale bar, 10 µm). (B) Reconstructed micropore structure. (Upper) Side view. (Lower) Top view. Color represents the absolute value of the radial parameter. Pores were located in 3D using an automated procedure, and areas of high slope (shown in gray) were filtered out. (C) Distribution of 3D nearest-neighbor distance. (D) A weak dependence of pore-to-pore distance on diatom size is observed (slope 0.003, significant at ). Each value is obtained by finding the peak of the distribution of interpore distance by fitting a Gaussian curve. (E) Top view of a reconstructed micropore structure, with color representing the phase of the radial parameter. (F) Histogram of the absolute value of . (G) Histogram of the phase of .

For a given valve, the spacing between pores was very regular. The distribution of nearest-neighbor pore spacings from a representative valve is shown in Fig. 1. In this example, the average spacing obtained is 1.4 µm.

While the overall size of valves varies by a factor of two (from 50 to 110 µm in diameter), the interpore spacing is very consistent, ranging from about 1.2 to 1.4 µm (Fig. 1). We observe a small but significant dependence of the pore spacing with the size of the valve (Fig. 1) (slope 0.003, significant at ). This suggests that pore size and spacing are predominantly determined by local processes (e.g., phase separation), mostly independent of those that control the overall valve dimensions.

Micropores show an intricate arrangement in the valve. In the vicinity of a single pore, the structure closely resembles a triangular lattice. However, a straight line drawn in any direction from the center to the edge of the valve will tend to follow the lattice.

To quantify this unusual structure, we adapted a geometrical tool developed for studies of two-dimensional (2D) crystals (e.g., ref. 25). The hexatic order parameter, , is defined with the complex expressionwhere l is a point on the lattice and m is one of its n nearest neighbors (determined by Delaunay triangulation). θ is the angle between the line connecting the lth and mth points and a reference axis. The amplitude of is 1.0 for a perfect triangular lattice, but for a random arrangement of points, we observe a distribution peaked at 0.3. Therefore, we can use to locally assess crystalline order around each pore. This is visualized for one valve in Fig. 1, where each pore is colored according to . A histogram of these values is shown in Fig. 1.

Most of the values for this valve are close to one, but there is a very broad tail of low values. These more disordered areas are more frequent near the center of the valve (Fig. 1).

On the other hand, the phase of the hexatic order parameter, , shows the local orientation of the lattice relative to a selected axis. Usually, a fixed axis is selected as a reference. Since our structure appears to be radially aligned, we found it more useful to select , the unit vector connecting the center of the valve to the pore, as the reference axis for each point. By construction, this unusual choice of reference axis has no effect on . However, it clearly shows how the local orientation of the crystal compares with . These orientations are shown point wise for a typical valve in Fig. 1. The distribution of , shown in Fig. 1, is sharply peaked at zero. This demonstrates that the nearest neighbors of each micropore tend to be arranged in a hexagon that is oriented with one vertex pointing toward the center of the valve.

Thus, the pattern of holes in the valve can be thought of as a lattice constantly turning around the center of the valve. We define the overall crystallinity of a valve, . This will equal one for a perfect crystal. We define the overall radial alignment of a valve, . A well-aligned valve has narrow distribution of angles, and therefore, a relatively large value of .

A scatterplot showing the crystallinity and alignment parameters of N = 86 valves is shown in Fig. 2. Note that the values of and are correlated. This is counterintuitive, as crystallinity and radial alignment are incompatible. A closer look at the data shows that more highly organized valves tend to be larger, as indicated by the color coding in Fig. 2.

Fig. 2.

C. granii valves display a range of alignment, crystallinity, and defect patterns. (A) Scatterplot of alignment () and crystallinity () derived from , with a gray oval to guide the eye. Data point color shows diatom size in micrometers. Larger diatoms tend to be both more aligned and more crystalline. (B) Example of a confocal stack projection of the diatom labeled in A. (Scale bar: 10µm.) (C) Reconstructed pores of the same diatom, showing points with six neighbors in gray, five neighbors in red, and seven neighbors in blue. Pairs of five to seven (red–blue) defects are scattered across the structure, with five (red) mostly located closer to the edge and seven (blue) toward the center. (D) Crystal lines traced over a section of a confocal frustule image. Defects are shown in red (five) and blue (seven). The radial direction is indicated in green. (E) Cumulative number of defects (sevenfold is shown in blue, and fivefold is in red) plotted against their radial position. The dashed line shows a slope of 7.3. (F) A histogram of slopes for all observed diatoms, with same color scheme. The vertical dashed line is showing . (G) Photograph of a crochet circle, with double stitches shown in white.

To further clarify the competition between alignment and crystallinity, we identified defect sites in the lattice. Using Voronoi tessellations, we determined the number of nearest neighbors for each pore. In a triangular lattice, each point has six nearest neighbors. In a real 2D crystal, there can be a significant number of defect sites with five or seven neighbors. The defect pattern of an exemplary valve is shown in Fig. 2 (defect patterns for additional valves are found in ). There, points with six neighbors are shown in light gray, while five- and sevenfold defects are highlighted in red and blue, respectively. Most defects are found in five–seven (red–blue) defect pairs. These pairs are mostly isolated, only occasionally forming higher-order structures. This is in contrast to 2D crystals near equilibrium, where defect pairs tend to condense into chains, forming boundaries between crystal grains.

Intriguingly, defect pairs are polarized, with the sevenfold defect closer to the center than the fivefold one. Note that five–seven defect pairs correspond to the termination of a crystal line, known as an edge dislocation, illustrated in Fig. 2 (26). To avoid overlapping pores, crystal lines must terminate as they converge on the center of symmetry. The fivefold defect sits at an end of the terminated crystal line, and the sevenfold defect is found just beyond it. The observed polarization of the defect pairs is, therefore, a direct consequence of radial alignment of the lattice vectors.

Defects are more frequent near the center of the valve. To quantify this, we plot the cumulative number of defects within a distance, r, from the center for an exemplary valve in Fig. 2. For a conventional 2D crystal, the defects are distributed homogeneously, and their cumulative number would increase as r2. However, for diatom valves, the cumulative number increases with linearly with r. The mean slopes are 6.8 ±0.1 and 6.3 ±0.1 defects per lattice spacing for five- and sevenfold defects, respectively (Fig. 2). Note that a linear scaling of the cumulative number of defects is indicative of a defect density scaling like .

While radially aligned crystal structures are, to our knowledge, unknown in materials, they are common in the traditional craft of crochet, as shown in Fig. 2. In a circular crochet, the square crochet lattice grows radially because each new stitch is inserted, by default, in a stitch on the previous row. However, every crocheter knows that each new row must be expanded by approximately six stitches for the circle to lay flat (27). This is geometrically necessary, as the circumference of each layer increases by for each row, where a is the spacing of the stitches.

The same argument applies to the diatom lattice, but the increase in the number of lattice sites per row is slightly different because the underlying lattice is triangular (). When the triangular lattice is aligned with the radial direction, the defect insertion rate is defects per lattice spacing. This value is comparable with the range of values we observe, as shown in Fig. 2.

So far, we have focused on structural analysis of a single species, C. granii. The morphology of diatom frustules, however, is tremendously diverse. Within the Coscinodiscus family, frustules are generally barrel shaped with flat round valves featuring organized patterns of microscopic pores. We collected 2D images of the valves of other Coscinodiscus diatoms from a variety of sources, shown in Fig. 3. We located pore positions using a semiautomated procedure (Fig. 3) and analyzed them the same way as our C. granii data. The crystallinity and alignment parameters for these diatoms are shown in Fig. 3. While samples from Coscinodiscus radiatus and Coscinodiscus wailesii fall into the range of morphological parameters we observed in C. granii, samples from Coscinodiscus angustelineatus and Coscinodiscus eccentricus are completely different. They show excellent hexagonal order, with no indication of radial alignment.

Fig. 3.

Structural diversity of Coscinodiscus family diatoms. (A) Images of diatoms. C. angustelineatus and C. eccentricus () are from ref. 28, Coscinodiscus curvatulus is from ref. 29, and C. wailesii and C. radiatus were provided by Johannes Goessling, Centre for Environmental and Marine Studies, University of Aveiro, Aveiro, Portugal. (Scale bars: 10µm.) Reprinted with permission from refs. 28 and 29. (B) Corresponding points obtained using the pore locating procedure. (C) Scatterplot of alignment () and crystallinity () derived from , with the gray oval representing the observed range of C. granii values. (D) Confocal stack average projection of a C. granii diatom grown in the presence of germanium. (Scale bar: 10µm.) (E) Corresponding points obtained using the pore-locating procedure.

For C. granii, alignment of pores is very robust, even when crystallinity is destroyed by elemental substitution. The structure of the silica frustule can be altered by growing diatoms in the presence of germanium (30, 31). Since germanium and silicon belong to the same group of the periodic table, they have similar chemistries. When germanium is introduced into the cell-culture medium, it is taken up by the cell and incorporated into the frustule. Over a range of germanium concentrations, diatoms continue to grow and divide at a normal rates (). However, the morphology of the frustule can change dramatically, as shown in Fig. 3 and . Laser ablation inductively coupled plasma time-of-flight mass spectrometry (LA-ICP-TOFMS) indicated insertion of germanium at a uniform rate across the frustule (). While the morphology of germanium-perturbed valves changes significantly from cell to cell, the pores tend to have a wider range of sizes and shapes, which include crack-like gaps in the frustule. In these conditions, crystalline order has essentially disappeared. Nevertheless, strong alignment of pores along the radial direction is evident in the micrographs (e.g., Fig. 3).

We have found that the organization of pores in diverse natural and perturbed valves is characterized by a blend of crystallinity and radial alignment. On one hand, natural variations across the Coscinodiscus family maintain crystallinity but feature varying degrees of radial alignment. On the other hand, perturbation of C. granii with germanium appears to maintain radial alignment while destroying crystallinity. This suggests that crystallinity and alignment are driven by different mechanisms.

Since the biochemical mechanisms that regulate these patterns are unknown, we developed a minimal physical model to investigate the competition of radial alignment and crystallinity during morphogenesis of the frustule. Following Sumper (7), we modeled the rearrangement of phase-separated domains that template the silica. Crystallinity naturally emerges due to repulsive interactions of the domains. To favor alignment, we introduced an ad hoc torque on interacting pairs of domains. The mechanism of alignment in the frustule is unknown. We suggest a few possibilities in the conclusions. Specifically, we considered 2D systems of N Brownian disks with radius σ, placed randomly within a circular confining potential. The disks repel each other with a force proportional to their overlap, quantified with a spring constant k. An alignment force tends to rotate the bond between neighboring particles to align with the local radial direction. The magnitude of the alignment force also depends on the degree of overlap and is proportional to a spring-like constant, k. More details are in .

We performed a series of Brownian dynamics simulations with a varying strength of alignment, k, and fixed repulsion, k. Particles were randomly distributed within the circular simulation domain with an effective area fraction and evolved in time until reaching steady state. Two system sizes were considered, with and , roughly matching the smallest and largest numbers of pores observed for the valves of C. granii. The resulting steady structures are shown in Fig. 4 , and additional structures can be found in .

Fig. 4.

A minimal numerical model yields diatom-like structures. (A–C) Pore locations from simulations. In Left, Center, and Right, k = 0 and and , and and , respectively. (A) Color represents the amplitude of . (B) Color represents the phase of . (C) Color represents the numbers of nearest neighbors. (D) Scatterplot of alignment () and crystallinity () derived from , with the gray oval representing C. granii data. Triangle color represents the value of k used in the simulation. For each value, the simulation was run five times, and for each run, the initial distribution of points was randomized. (E) The cumulative number of defects plotted against the distance from center for k = 0 and 0.8 with . Fivefold defects are plotted in red, and sevenfold defects are in blue. The dashed line shows a slope of 7.3.

The crystallinity and alignment parameters of the resulting structures are shown in Fig. 4. As k increases, the crystallinity parameter drops, and the alignment parameter increases. While alignment increases only modestly for the smaller system size (triangles), it increases much more strongly for the larger system size (circles). At large values of k, the simulated systems achieve values of crystallinity and alignment typical of C. granii valves. For a fixed value of k, the larger valves show better crystallinity and alignment parameters, echoing the size dependence of real valves, shown in Fig. 2. The cumulative counts of defects as a function of the distance from the center for k = 0 and 0.8 are shown in Fig. 4. Like in real diatoms, the simulated structures at accumulate defects linearly with increasing distance from the center. Here, the slope is very close to the expected value of , indicated by the dashed line in Fig. 4. Defect pairs in the simulation tend to organize into higher-order structures separating regions with homogeneous alignment. These grain boundaries are not observed in the real C. granii valves, suggesting that the pore structure may have kinetically arrested before finding its energetically preferred configuration (5). Earlier time points in the simulation more closely resemble the structure of the valve, with similar levels of crystallinity and alignment, but no clear grain boundaries (Movie S1).

Conclusion

We have found that micropores of C. granii diatom valves are arranged on radially aligned triangular lattices. These beautiful structures feature a geometrically necessary increase of defect density near the center of the valve. While triangular lattices are observed across the Coscinodiscus family, the degree of radial alignment varies from species to species. When the crystallinity of C. granii is perturbed by the introduction of germanium, the valve structure remains radially aligned. These observations suggest that these two structural features have independent physical origins. The range of structures observed across the Coscinodiscus family can be recapitulated with a minimal model that includes conventional steric repulsions as well as a torque aligning their separation vectors with the local radial direction.

While steric repulsions are consistent with the self-assembly of stable phase-separated domains, the origin of the torques producing alignment is unknown. At a minimum, this requires the imposition of a center of symmetry, which may be achieved by actin and microtubules, whose filaments have been shown to be radially aligned in close association with the outer surface of the SDV (20). Alternatively, radial alignment could be imposed from the edge of the circular domain (32) or by organization along radially aligned creases of the SDV due to confinement within the curved boundary of the cell (33, 34). Direct live-cell observations of frustule development will be essential to identify the mechanisms driving alignment and more generally, to understand the phase separation process that produces micropores.

Whatever its origin, radial alignment imposes geometrically necessary defects in a flat 2D lattice. Specifically, pairs of five- and sevenfold defects need to be inserted with a frequency inversely proportional to the distance from the center of symmetry. This is in stark contrast with 2D crystals conforming to a curved surface in the absence of a center of symmetry. There, five- and sevenfold defects are introduced in an imbalanced manner, with a topological charge density proportional to the local Gaussian curvature (35, 36).

Traditional knowledge of crochet demonstrates that radial alignment in the absence of defects imposes curvature (37). Indeed, a defect-free radially symmetric crochet lattice deforms to produce a tube. Furthermore, an appropriate insertion of defects produces fascinating structures, including hyperbolic planes (38) and Lorentz surfaces (39). More generally, insertion of defects could provide a route to control the embedding of a 2D object in 3D space by specifying the local curvature (40–42).

Diatoms are master builders of inorganic structures at the microscale. We expect that a deeper understanding of the processes they use to sculpt biomineralization could inspire approaches to sustainable micro- and nanomanufacturing.

Materials and Methods

C. granii Cells.

C. granii (strain K-1834) diatoms were purchased from The Norwegian Culture Collection of Algae (Oslo, Norway) and cultured in L1 media (43). The cells were cultured in filter-top flasks in a temperature- and light-controlled incubator (Tritech Research, Inc.) at 15 ∘C under a red and blue light-emitting diode grow light. The light was cycled 14/8 h on/off, and light flux when on was 200 to .

Fluorescent cell walls were obtained by metabolic insertion of Rhodamine B (Sigma-Aldrich). The Rhodamine dye was added to the culture medium at a concentration of 1µM, and the culture then was grown for at least 24 h (44). Diatom cell walls were extracted using a “soft-cleaning” method based on ref. 44. Cells were collected by centrifugation at 2,400 × g, resuspended in 0.1M Na–ethylenediaminetetraacetic acid (EDTA) and 2% sodium dodecyl sulfate solution, and heated at C for 5 min.

Imaging.

Imaging was performed with a Nikon eclips-Ti microscope (Nikon Corporation) equipped with a spinning disk confocal attachment (Intelligent Imaging Innovations, Inc.). Imaging was performed with a Nikon PLANAPO 100× NA 1.45 lambda oil immersion objective or Nikon PLANAPO 100× NA 1.45 lambda oil immersion objective (Nikon Corporation) for germanium diatoms.

Image Analysis.

Image analysis of C. granii diatoms was performed in MATLAB (MathWorks) using an automated procedure. The image z stacks were contrast adjusted, inverted, preprocessed with a band-pass filter, and binarized using the MATLAB function imbinarize. For each object, area and solidity were determined using regionprops, and both were used for filtering. Z positions were determined by first finding the track of each object through the z stack following the method in ref. 45 implemented for MATLAB in ref. 46. Then, for each object, the z plane with the largest area was selected as the final z coordinate. The center of the structure was determined manually.

Image analysis of other Coscinodiscus species was performed on 2D images using MATLAB. The images were preprocessed by being inverted (when necessary), filtered with a band-pass filter, thresholded, and binarized. The function regionprops was used to locate features and remove noise by thresholding feature size. The center of the structure was determined manually.

Elemental Analysis.

LA-ICP-TOFMS element imaging (47) was performed using the single-pulse mode according to Neff et al. (48). A laser ablation system (GeoLas C; Lambda Physik) employing a modified tube cell design (49, 50) was coupled to an inductively coupled plasma time-of-flight mass spectrometry instrument (icpTOF2R; Tofwerk AG). The piezo–electrically driven xyz translational stage (1-nm scan resolution; SmarAct GmbH), the laser, and the time-of-flight mass spectrometry data acquisition were triggered using a custom-built imaging control system (48). The full mass spectrum from Na to U was acquired, and the quantification of element mass fractions was performed for each pixel based on 100% mass normalization (51) using National Institute of Standards and Technology standard reference material (NIST SRM 610) as an external reference material and assuming that all elements are present as oxides.

Numerical Model.

Consider a set of N disks with radius σ, located at , where . The equation of motion of a single object is given bywhere b is the hydrodynamic mobility.

Disks interact with a thermal bath, which generates a random force, . This force is uncorrelated in space and time and drives Brownian motion, characterized by a diffusion coefficient D. The disks are confined to a circular domain with radius R by an external force

The center of this circular domain is taken to be the origin for all of the particle locations.

Overlapping disks interact in two ways:where is the force exerted on the ith particle due to the jth particle and is the vector separating particles i and j. The first contribution to the interaction is a conventional spring-like repulsion,

whose strength is determined by a spring constant, k. With only this interaction term, sufficiently concentrated particles will evolve into a crystalline steady state with a triangular lattice. Therefore, we introduced a second interaction that works to point the separation vector toward or away from the center of the structure,

where This force is applied perpendicular to the separation vector along , where is the unit vector pointing out of the plane. Like the repulsive force, it depends on a spring constant, k, and is proportional to the extent of overlap between the particles. However, instead of pushing the particles away from each other, it rotates them around their center of mass, . The strength of increases as turns away from . This specific form is not motivated by a specific physical model but is simply meant to capture a generic tendency for alignment.

Simulations were performed in MATLAB (MathWorks). We give the disks σ = 1, b = 1, and D = 0.01. The filling fraction was kept constant at . Simulations were run with a time step of for a total number of time steps . The radius of the confining wall was R = 40 or 85 for or particles, respectively. The stiffness of the particles and wall was fixed at . The alignment parameter k was varied from 0 to 1.2. The simulations were run five times for each value of k, with initial particle positions randomized for each run.