Significance Rangeomorph fronds characterize the late Ediacaran Period (575–541 Ma), representing some of the earliest large organisms. As such, they offer key insights into the early evolution of multicellular eukaryotes. However, their extraordinary branching morphology differs from all other organisms and has proved highly enigmatic. Here we provide a unified mathematical model of rangeomorph branching, allowing us to reconstruct 3D morphologies of 11 taxa and measure their functional properties. This reveals an adaptive radiation of fractal morphologies which maximized body surface area, consistent with diffusive nutrient uptake (osmotrophy). Rangeomorphs were adaptively optimal for the low-competition, high-nutrient conditions of Ediacaran oceans. With the Cambrian explosion in animal diversity (from 541 Ma), fundamental changes in ecological and geochemical conditions led to their extinction.

Abstract The branching morphology of Ediacaran rangeomorph fronds has no exact counterpart in other complex macroorganisms. As such, these fossils pose major questions as to growth patterns, functional morphology, modes of feeding, and adaptive optimality. Here, using parametric Lindenmayer systems, a formal model of rangeomorph morphologies reveals a fractal body plan characterized by self-similar, axial, apical, alternate branching. Consequent morphological reconstruction for 11 taxa demonstrates an adaptive radiation based on 3D space-filling strategies. The fractal body plan of rangeomorphs is shown to maximize surface area, consistent with diffusive nutrient uptake from the water column (osmotrophy). The enigmas of rangeomorph morphology, evolution, and extinction are resolved by the realization that they were adaptively optimized for unique ecological and geochemical conditions in the late Proterozoic. Changes in ocean conditions associated with the Cambrian explosion sealed their fate.

In parallel with large-scale geochemical transitions associated with ocean oxygenation (1⇓–3), the Ediacaran Period (635–541 Ma) records a major diversification of multicellular eukaryotes. Rangeomorph fronds (575–541 Ma) dominated early Ediacaran biotas (4) and have a characteristic branching morphology, distinct from any known Phanerozoic organism (5). Although the fronds are often preserved as flattened impressions, exceptional moldic fossils preserve details of the 3D branching structure to a resolution of 30 μm (6). Qualitative classifications for rangeomorph branching patterns have been proposed (7, 8), but no quantitative model has previously been formulated. Because branching is repeated over decreasing size scales (with up to four observed orders of branching), rangeomorph fronds have been informally described as self-similar and fractal (4⇓–6). Although this has potential implications for the functional optimality of their morphologies (5, 9), the extent to which they are formally fractal and self-similar (10, 11) has not previously been tested. Furthermore, until now, evolutionary transitions in branching patterns have not been characterized within any quantitative framework.

Rangeomorphs inhabited shallow to abyssal marine environments (1, 8, 12⇓–14), evidently precluding photosynthesis for most taxa (12). Preservational features, including bending and overfolding (4, 15), suggest that rangeomorphs were soft-bodied. No evidence exists for either motility or active feeding (such as musculature, filter feeding organs, or a mouth). Consequently, rangeomorphs have been reconstructed as sessile, feeding on organic carbon by diffusion (or possibly endocytosis) through the body surface (3, 5, 12, 16), with a large surface area to volume ratio aiding nutrient uptake (5, 16). The adaptive potential of their different branching morphologies has, however, never been quantified.

Here, using parametric Lindenmayer systems (L-systems) (17, 18), we present a unified model to describe the branching structure of Ediacaran rangeomorph fronds. Our quantitative parameters provide a far more detailed definition of frond morphology than was previously possible. These parameters are then used to reconstruct 3D space-filling strategies within the group and evaluate potential frond functions, revealing an adaptive radiation of fractal organizations. This provides a new framework for the study of growth, functional morphology, and evolution of these lost constructions.

Conclusions With their terminal Proterozoic extinction and unique morphology, Ediacaran rangeomorph fronds have been described as a failed experiment in evolution (24). However, our analysis demonstrates that these intriguing fossils possessed a fractal morphology which combined programmatic [and potentially genetic developmental (9, 11, 28)] simplicity, structural versatility, and functional optimality for the uptake of organic carbon (osmotrophy) and oxygen. The appearance of rangeomorph fossils (1) occurred after a move away from anoxic, sulfidic, and ferruginous oceans, toward conditions more favorable for aerobic macroorganisms (2, 29). Their disappearance coincides with the Cambrian explosion in metazoan diversity, a dramatic increase in competition, and, crucially, decreased availability of organic carbon in ocean water (2, 12, 19, 24, 34, 36). These potentially interacting factors suggest that the Ediacaran to Cambrian transition was a bad time to be a sessile, soft-bodied osmotroph. The unique rangeomorph fronds were fractal, surface area specialists of the Ediacaran. At the Cambrian explosion, the ecological and geochemical conditions to which the rangeomorphs were optimized ceased to exist, and their extraordinary body plan was lost from life’s repertoire.

Methods First, this study established a unified model for rangeomorph theoretical morphology using parametric (18) Lindenmayer (L) systems (17), written within the L-studio programming environment (37). L-systems are a class of parallel derivation grammar, in which specified production rules are applied in parallel to control iterative rewriting of the axiom (a starting string, here representing the first stem segments, and holdfast if present). The symbols produced are then interpreted graphically to visualize a geometry encoded by the output L-system string (representing the branching system) (18). A branching L-system is characteristically fractal, with self-similar elements visible at decreasing size scales (38). Parametric L-systems allow branching parameter values (such as branching angles and growth rates) to vary between branches (for example, of different orders or ages), enabling realistic representation of biological structures with nonuniform branching patterns (18). Rangeomorph morphologies were modeled using formal production rules and parameters based on branching patterns of 11 studied species and quantitative measurements (of branching angles and body dimensions) from best-preserved representatives (SI Appendix). For each species, the L-system output geometry was then processed and analyzed using Blender 2.69 and MATLAB (2012b; Mathworks). For measurement comparability, L-system output was standardized to the same finite approximation (four orders of branching). The fractal dimension of each modeled frond was estimated using box counting, the method most commonly used to analyze complex fractal shapes (10). This method determines the number (n) of boxes of a given size (r) that cover the input image. The scale (r) is incrementally decreased to determine the relationship between box size and image coverage. The slope of the line for this relationship gives the fractal dimension of the image D = −log(n)/log(r). If dimension D is not an integer, this indicates that the image is a fractal (its geometry does not correspond exactly to an integer dimension, i.e., a 1D line, 2D plane, or 3D volume). The input image for 2D box counting was a 2D binary (black or white) skeleton of the branching pattern, in frontal view (SI Appendix, Fig. S13). Input for 3D box counting was a 3D binary skeleton of the branching pattern (SI Appendix, Fig. S13). Two-dimensional box counting was conducted using ImageJ (39). Three-dimensional box counting was conducted in MATLAB with a script incorporating the Wavefront object toolbox (40), Inhull function (41), and boxcount toolbox (42). Because 3D box counting is highly computationally intensive, two large species, Fractofusus misrai and Pectinifrons abysallis, were analyzed using 2D box counting only. Functionally relevant frond properties were calculated from the output mesh, including the surface area to tissue volume ratio [with modeled tissue depths of 0.1, 0.5, and 1 mm (comparable to ref. 5)] and the size (height, width, and depth) of a bounding box around the modeled organism. Bounding box axes for each species were oriented relative to inferred life position (based on fossil morphology and preservation). Dimensions for each species were scaled based on approximate specimen lengths recorded in the literature (SI Appendix, Table S2). A Euclidean cluster analysis of these scaled dimensions was conducted using paired group linkage in palaeontological statistics (43).

Acknowledgments We thank A. Liu, J. Gehling, C. Kenchington, and M. D. Brasier for discussion and making fossil casts and photographs available for study, with additional thanks to the Oxford University Museum of Natural History and the University of Cambridge Sedgwick Museum of Earth Sciences for access to their collections and to the South Australia Museum for providing specimen photographs. We also thank S. Gerber and three anonymous reviewers for highly constructive comments on the manuscript. This research was supported by Templeton World Charity Foundation Grant LBAG/143.