Standardizing Estimates of Biomass at Recruitment and Productivity for Fin- and Shellfish in Coastal Habitats

Assessing the biomass and productivity of fin- and shellfish supported by coastal ecosystems is important to develop plans for the conservation and restoration of these ecosystems, but such assessments are not easy to obtain. We developed a protocol that, from density data, quantifies biomass at recruitment for species where information exists to derive life history tables, and productivity where such information does not exist. Our protocol also assesses the variability (i.e., variance) for the calculated biomass and productivity values. For relatively well-reported species, inferences regarding differences among habitats or species can be suggested. For instance, application of our protocol to juvenile pinfish confirms its well-known preference for structured habitats. Mud crabs also seem to reach higher productivity levels in structured than open bottom habitats. For poorly reported species, only a general idea can be gleaned. However, larger data sets of fin- and shellfish density in shallow coastal systems are needed to increase the accuracy, precision, and comprehensiveness of the estimates of biomass at recruitment and productivity generated with our protocol. With such larger data sets and the use of statistical tools such as Bayesian methods, the protocol can significantly help improve our understanding and management of fisheries productivity in coastal systems.


Introduction
Coastal ecosystems lie at the interface between the land and ocean, and include estuaries and other types of coasts that are not at the confluence of rivers and oceans (Valiela 2006). Foundational components, such as marshlands, seagrass beds, biogenic reefs, and sediment flats, form these ecosystems. Coastal ecosystems generate many benefits for humankind and wildlife. They provide habitat for a plethora of species, including fin-and shellfish, birds, and mammals. They can also act as filters of land-derived pollution before it enters the open ocean. Primary producers, through nutrient uptake, and microorganisms, through denitrification and other microbial processes, may significantly reduce nutrient loading from land into the open ocean (Tobias et al. 2001, Sparks et al. 2015. Coastal ecosystems can also be large carbon sinks and help mitigate anthropogenic CO 2 accumulation and climate change (Duarte et al. 2005, Fourqurean et al. 2012. Coastal ecosystems can buffer wave energy substantially and protect the coastline against storms (Manis et al. 2015, Sharma et al. 2016. The habitat role of coastal ecosystems has received much attention. Many species of fin-and shellfish with commercial and recreational importance utilize these ecosystems. Most of such species, however, reside in the shallower parts of coastal ecosystems as juveniles and pre-adults, while adults move to deeper parts of the coastal ecosystem, or even to farther offshore waters, and do not return to the shallow waters to which they recruited ). In some other species, adults may temporarily (e.g., a few months per year) return to shallow waters (Sheaves et al. 2015). More than 90% of all commercial fish landings in the Gulf of Mexico are species that use shallow coastal waters during some portion of their life cycle (Kennish 1999, Lellis-Dibble et al. 2008). Exported production from shallow to deeper waters through juvenile migration represents prey for fish in the deeper waters. In turn, the invertebrate and fish species that reside permanently in the shallow waters of coastal ecosystems constitute prey for juveniles, pre-adults, and adults of temporary residents. All together, these various species of fin-and shellfish that occur in shallow waters of coastal ecosystems represent a prolific trophic resource for large apex predators that frequently visit these waters, such as sharks, and for permanent and migratory waterfowl, many of which are targeted for conservation and recreation (e.g., bird watching). Endangered megafauna, such as manatees and sea turtles, may also forage in the shallow parts of coastal ecosystems (Aven et al. 2015). Therefore, coastal ecosystems have many economic and societal ramifications and are often the object of highly debated management policies.
The need for effective management of coastal ecosystems has been heightened by the major losses and degradation that these ecosystems have endured in recent decades. Human activities are often the cause of such degradation. As humans develop the coast, maritime forests and wetlands are removed and replaced with urban landscapes. The conversion of natural into developed coastal land often results in increased sediment and nutrient inputs into coastal waters, which can severely impair water quality and cause the decline of important habitats such as seagrass beds and biogenic reefs (Valiela 2006). Human overexploitation of coastal resources may accelerate the degradation of these habitats. Past reports have estimated coastal habitat losses amounting to 29% for seagrass beds (Waycott et al. 2009), 67% for marshlands (Lotze et al. 2006), and 85% for oyster reefs (Beck et al. 2011) over large areas across the world. Large losses of these habitats have also been reported for the Gulf of Mexico (Handley et al. 2007, Zu Ermgassen et al. 2012, Sparks et al. 2013. Much effort is ongoing to offset these losses through conservation, restoration, and mitigation practices.
Two key metrics to gauge the habitat value of coastal ecosystems are fin-and shellfish biomass (generically defined as mass per bottom area unit, for instance dry weight per square meter of bottom) and productivity (generically defined as the gain in biomass over time, for instance dry weight of fish biomass produced in a year per square meter of bottom) in these ecosystems. This information can help develop management policies for coastal ecosystems, for instance by setting quantitative estimates of fisheries productivity that can be maintained or enhanced through habitat conservation and restoration (Lellis-Dibble et al. 2008, Sheaves et al. 2015. However, these metrics are difficult to measure, their availability in the published literature is limited, and the information available covers varying and often discrepant taxonomic diversity, habitats, time periods, and gear types. In addition, many available estimates are subject to several assumptions and, in some instances, considerable uncertainty that is often not well quantified (Peterson et al. 2003, French McCay et al. 2015. This limits our understanding of the habitat role of coastal ecosystems as well as our ability to efficiently manage these ecosystems.
In this paper, we develop a protocol that, departing from density values (expressed in number of individuals per square meter of bottom), derives two kinds of estimates. The first kind is estimates of biomass at recruitment. In these calculations, recruitment corresponds to the post-larvae or early juveniles that settle in the shallow parts of coastal ecosystems. The calculations can be applied to permanent resident (i.e., all life cycle stages-eggs, larvae, juveniles, and adults remain in the shallow parts) and temporary resident species (i.e., large juveniles, pre-adults, or adults migrate to deeper coastal or farther offshore waters) of shallow parts of coastal ecosystems. Biomass at recruitment (expressed in grams of dry weight per square meter of bottom) corresponds to the biomass of post-larvae or juveniles of the given species that settle in the shallow coastal system at the time of recruitment. The calculation of estimates of biomass at recruitment following our protocol is contingent upon information provided by life history tables, which can be derived for most finfish species, and for some shellfish species, that permanently or temporarily utilize the shallow parts of coastal systems (Jensen et al. 1988, French McCay et al. 2003a, French McCay et al. 2015. However, for a few fin-and many shellfish species that permanently or temporarily utilize the shallow parts of coastal ecosystems, information to initialize and validate life history tables, and thus derive them with rigor, does not exist. In these cases, the protocol generates estimates of productivity for the given fin-or shellfish species. These estimates (expressed in grams of dry weight per square meter of bottom per year) correspond to the new biomass generated by the species over one year per bottom square meter in the shallow coastal system. The productivity estimates include the new biomass generated by all different species' life cycle stages as they naturally occur in the system. For both types of estimates (i.e., biomass at recruitment and productivity), we borrow from established methods to build a cohesive construct for their derivation, and importantly we also present how to quantify the uncertainty of the estimates. Our protocol allows for the derivation of estimates of biomass and productivity of fin-and shellfish species targeted for economic, recreation, or conservation purposes in coastal habitats, as well the quantification of uncertainty around such estimates. It can also help develop estimates of fisheries productivity under various scenarios of coastal habitat conservation and restoration. Thus, the protocol presented here can improve our capacity to effectively manage coastal ecosystems.

Derivation of Life History Tables
Here, we focus on the first year of life for the estimation of biomass at recruitment. Our protocol relies on comprehensive and detailed information of how fish size and mortality evolve throughout that first year of life (i.e., daily or, at most, weekly intervals). However, comprehensive records of speciesspecific measurements of growth and mortality rates during the first year of life directly obtained in the field are seldom available. Instead, where adequate information exists for model initialization and validation, detailed life history tables for the first year of life can be derived using well established equations. Life history tables assess how fin-and shellfish individual length, mass, and mortality rates change with individual age. The equations are different for larval and juvenile stages, and they are combined to provide estimates of daily growth (length) and mortality rates throughout the first year. Standard weight-length conversions are used as appropriate.
Here, we utilize estimates of mortality rates (as per day) and individual fresh weight (as grams of fresh weight per individual) over daily time intervals throughout the first year of life. The derivation of these estimates is explained in Appendix 1 (for further elaboration see French McCay et al. 2015) and the actual estimates provided in Appendix 2.

Species with Derived Life History Tables
We have chosen two fish species as examples, pinfish (Lagodon rhomboides) and black drum (Pogonias cromis). Species-specific values for model initialization and validation to derive the life history tables are obtained from Nelson (2002) for pinfish and Murphy and Muller (1995) for black drum. These species represent fish that recruit to shallow parts of coastal ecosystems, where they reside through their pre-adult stage and subsequently migrate to deeper waters as they become adults. Numerous records of fish density in coastal waters exist for pinfish; however, the number of density records for black drum in coastal waters is drastically lower. Thus, this comparison allows us to gauge how our calculations fare for well vs. poorly reported density data sets.
Pinfish is a widespread species in coastal ecosystems extending from Massachusetts (although rare north of Maryland) to the Yucatan Peninsula, and to Bermuda, and the Gulf of Mexico (Muncy 1984). Juveniles typically reside in the shallow parts of coastal systems from late winter to late fall, preadults move to deeper waters of coastal systems, and adults move farther offshore. The species is not commercially harvested, but it serves as an important consumer of invertebrates (as juveniles) and plants (as large juveniles and adults; Hoss 1974, Stoner 1982, and as prey (both as juveniles and adults) for harvested fish species (Jordan et al. 1996, Nelson et al. 2013). Pinfish has been extensively studied and its density in coastal ecosystems along the Gulf of Mexico is well documented.
Black drum inhabits a wider geographical range than pinfish extending from the Bay of Fundy to the North Atlantic and Gulf coasts to the South Atlantic coast (Argentina) (Sutter et al. 1986). The species occur in coastal habitats along this range, and it is common in the Gulf of Mexico. Juveniles stay in the shallower parts of coastal ecosystems from mid spring to mid fall, pre-adults (typically up to 2 years old) reside in deeper areas of coastal waters, and adults may move to even deeper areas farther from the coastline (Osburn andMatlock 1984, Cody et al. 1985). The species plays significant trophic roles as prey and predator. In addition, it constitutes important commercial and recreational fisheries in the Gulf of Mexico (Leard et al. 1993). In spite of this, there are relatively few reports of density for this species in coastal ecosystems along the Gulf of Mexico.

Species Without Derived Life History Tables
We have chosen two invertebrate species as examples, the mud crab (Rhithropanopeus harrisii) and the Gulf stone crab (Menippe adina). Life history tables cannot be derived with rigor for these species since information to initialize and validate them is not available. These two species reside in shallow parts of coastal ecosystems throughout their entire life cycles. Density reports for the Gulf of Mexico are more numerous for the former species; thus, this comparison allows us to gauge how this second method fares with varying levels of density data availability.
Mud crabs are omnivorous scavengers and mostly feed on algae; small invertebrates such as amphipods, copepods, polychaetes and bivalves; seagrass detritus; and other dead organic matter. The frequency at which they feed and the quality of what they eat depend on the habitat and their diurnal cycle of activity and foraging (Hegele-Drywa and Normant 2009;Williams 1984). Mud crabs can be found in coastal environments throughout the northern hemisphere, and they are considered global invaders introduced through ballast waters and commercial oyster shipments.
Stone crabs occur on sediment bottoms, oyster reefs, and rock jetties in coastal ecosystems. Adults burrow in mud or sand while juveniles hide among rocks. Stone crabs are high-level predators in waters in the South Atlantic Bight, Caribbean (Western Atlantic stone crab, M. mercenaria), and northern and western Gulf of Mexico (Gulf stone crab, M. adina) (Williams and Felder 1986). Stone crabs are commercially fished in the southeastern United States and managed as one species (Gerhart and Bert 2008).

Density Data Set
The data used in this paper is part of an extensive data set of nekton abundance in shallow habitats of coastal ecosystems extending from Laguna Madre in southern Texas to the Caloosahatchee River in southern Florida presented in Hollweg et al. (2019). The compilation of the data set, including the databases searched, identity of the variables compiled, criteria applied for data selection, and how the data were extracted or calculated, is explained in detail in Hollweg et al. (2019). This is a companion paper in the Estuaries and Coasts special issue "Quantifying the Benefits of Estuarine Habitat Restoration in the Gulf of Mexico" organized by M. V. Carle and K. Benson. The data set contains mean density values, expressed in number of individuals per square meter of bottom, obtained for specific habitats and time periods as reported by the studies compiled. Here, we used density values compiled for pinfish, black drum, mud crab, and Gulf stone crab. Due to the shallow nature of the habitats included in the compilation and the life histories of these species, density estimates correspond to young of the year (YOY) for pinfish and black drum, and are inclusive of all ages for mud crabs and Gulf stone crabs.
At this point, it is important to emphasize that our protocol, both for the estimation of biomass at recruitment and productivity, can be applied to both temporary and permanent resident fin-and shellfish species in shallow coastal systems. For the estimation of biomass at recruitment, we focus on the first year of life, and thus, the density values used for this estimation must only represent YOY. This should be mostly the case for density data obtained for temporary resident species in shallow coastal systems, since most of the individuals of these species that occur in such shallow systems are YOY (such as the two examples used here, pinfish, and black drum). In contrast, density values obtained for permanent resident species in shallow coastal systems should include more life stages other than YOY. Thus, estimations of YOY density from the wider population density values obtained for permanent resident species in shallow coastal systems must be first carried out before deriving estimates of biomass at recruitment for these species using our protocol. Total population density values must be used in the estimation of productivity with our protocol, since those estimates correspond to the new biomass generated by all different species' life cycle stages as they naturally occur in the system. For permanent resident species in shallow coastal systems, density values should include most life stages and such values can be used in the derivation of productivity estimates using our protocol. This is, however, not the case for density values of temporary resident species in shallow coastal waters, where efforts to estimate the fraction of life stages missing, and thereby produce density estimates that include all life stages within and out of the shallow coastal systems, are needed before deriving productivity estimates using our protocol.
The compiled density data for the four species targeted here (pinfish, black drum, mud crab and Gulf stone crab) included the following habitats: "near" submerged non-vegetated areas (within 5 m from fringing shoreline), "far" submerged nonvegetated areas (farther than 5 m from fringing shoreline), submerged aquatic vegetation (SAV), oyster reefs, and marshes. We divided the non-vegetated sites between "near" and "far" to account for shoreline edge effects Turner 1994, Minello andRozas 2002). The primary intent of this paper is to demonstrate how our protocol can derive estimates of biomass at recruitment and productivity for fin-and shellfish species in coastal habitats. Additionally, we also suggest potential uses of these estimates such as comparisons across species and habitats with the ultimate goal of informing management decisions. The main purpose of such comparisons is to offer some illustrative examples of uses of our protocol, and thus, we have restricted the comparisons to natural habitats (i.e., habitats that were not ostensibly degraded by human activities and/or that had not been restored by humans) to keep the comparisons simple and consistent. Using our protocol for comparisons between natural and restored coastal systems is definitely a promising venue of work that should be explored in future efforts.

Density Meta-analysis and Corrections
We followed the meta-analytic approach presented in Hollweg et al. (2019). We summarize the steps of this approach and we refer the reader to Hollweg et al. (2019) for further consultation. First, following an imputation method, we estimated the standard error (SE) for the mean density values in the data set where it was not reported or we could not calculate it based on the information available in the paper. Briefly, we used the expected relationship between the sample mean and sample standard deviation (SD) to impute missing SE (Hilbe 2014). Sample SD was regressed against sample mean from the records compiled, and tests were conducted to ensure the regression obtained was robust (Quinn and Keough 2002). SD was estimated from the regression for records with sample mean but not SE provided. If the sample size was not reported, we set it to n = 1.
Second, we calculated a weighted average and associated SE for all density entries for the same species corresponding to the same combination of habitat, sampling time, and gear type using a fixed effect model: where mean i is the i th mean of a given combination of habitat, sampling time, and gear type, w i is the weight of the i th mean, SE i is the standard error of the i th mean, and SE density weighted average is the standard error of the density weighted average. In the calculation of w i , we did not include a random error term encompassing variability due to author bias (i.e., several entries generated by the same authors) or the inclusion of different populations for the same combination of habitat, sampling time, and gear type because we did not have a sufficiently large sample size in the four species targeted to test for such random effects with rigor (Hollweg et al. 2019). Third, we applied correction factors for gear selectivity, capture efficiency, and recovery efficiency. Gear selectivity corresponds to the range of individual fish sizes that can be collected by the gear given its characteristics. Some gears do not normally capture fish smaller than a minimum or larger than a maximum size threshold Rozas 2002, Baker andMinello 2011). Capture efficiency corresponds to the fraction of size-apt fish within the sampled area that are actually enclosed and captured by the gear. Indeed some fish that are within the catchable size range out-swim and escape the gear as it is being operated (Rozas and Minello 1997). Recovery efficiency corresponds to the fraction of captured fish that is actually recovered from the gear and processed. Not all captured fish are necessarily recovered, particularly in gears with a secondary removal method (Rozas and Minello 1997).
Given the individual size ranges included in the data set for the four species considered, and the minimum and maximum size thresholds of the gear types (enclosure, towed and passive) in the data set, selectivity corrections were only deemed necessary for black drum collections with enclosure-type gears. This is further elaborated in the "Results" section as we address each of the four species separately. In contrast, corrections for gear capture and recovery efficiency apply to all four species and gear types considered. Along with gear type, capture and recovery efficiency also depend on the habitat considered. Thus, we developed correction factors for capture and recovery efficiency for all combinations of habitats and gear types compiled in our data set. To do that, we carried out an extensive literature search and, for each combination of habitat and gear type in our data set, we derived a mean conversion factor and SE. We did this separately for capture and recovery efficiency. The procedure is detailed in Hollweg et al. (2019).
We then calculated overall gear efficiency for each combination of habitat and gear type as: where C hg and R hg are the capture and recovery efficiency for the h th habitat and g th gear type. The variance of G hg was calculated using the equation reported by Goodman (1960) that provides an unbiased estimate of the exact formula of the variance of the product of two independent random variables: The entire populations of the capture and recovery efficiency values for each combination of habitat and gear type are not known in their entirety and without uncertainty. Our efforts, as exhaustive as they may be, can only provide a number of values out of the entire populations of those values. Thus, the mean for capture and recovery efficiency and their variances are based on a limited sample and not on the entire population of values. Because of this, population moments need to be replaced by the corresponding sample moments, and the exact equation of the variance of the product of two independent random variables is converted into its unbiased estimate depicted in Eq. 5 (for further elaboration see Goodman 1960). The inevitably limited sample size in our calculations, as it is the often the case in ecological studies, implies that using the exact equation of the variance of the product of two independent random variables is not as accurate as using its unbiased approximation.
Overall gear efficiency for each combination of habitat and gear type was used to correct the density weighted averages for specific combinations of habitat, sampling time, and gear type: where using the unbiased estimate of the exact formula of the variance of the product of two independent random variables reported by Goodman (1960): Subsequently, we averaged all the density values corrected for overall gear efficiency that corresponded to the same combination of habitat and sampling time (D ht ): and calculated its variance as: where D G ht are the density values corrected for overall gear efficiency corresponding to the specific combination of habitat and sampling time, and N ht is the count (number) of such values. All these steps are common to our calculations of biomass at recruitment and productivity ( Fig. 1).

Calculations of Biomass at Recruitment
The next step for these calculations was, using life history tables, to estimate density at recruitment from values of mean density at a given sampling time post-recruitment (Fig. 1). The fraction of YOY present at the beginning of day 1 that remain alive at the end of the day (YOY 1 ) corresponds to: where m 1 is the mortality rate for day 1 and is expressed in day −1 . In turn, the fraction of YOY remaining after t days (YOY t ) corresponds to: where m i is the mortality rate for day i, m t is the mortality rate for day t, and YOY t − 1 is the fraction of YOY present at the beginning of day 1 that remain alive at the end of day t-1. All mortality rates are expressed in day −1 and can be obtained as modeled values provided in the life history tables (see Appendix 2 for actual values). Following this, the density at recruitment (D R ) can be estimated from the density obtained at a sampling time postrecruitment (D ht ) as: where t covers the time span elapsed from recruitment to sampling. This procedure allowed us to derive a separate estimate of density at recruitment for each density value in the data set obtained at a later date as specified by the sampling time.
To estimate the variance of D R , we followed a multi-step process. The first step was to derive variance estimates for the daily mortality rates. To do this, we used the linear regression model provided by Bradford (1992): where m d is the mean interannual daily mortality for day d, var (m d ) is the variance of the daily mortality values that compose the interannual mean, and ln denotes natural logarithm. The model was generated using a literature survey of mortality rates for egg, juvenile, and adult stages of marine, freshwater, and anadromous fish species. At least two values of daily mortality rates corresponding to different years were obtained for each species. The mean interannual daily mortality was calculated for each species, and after transformation to natural logarithms, the variance of the interannual daily mortality values was regressed against the mean (see Bradford 1992 for further details). Thus, this effort includes temporal variability in the estimates of variance for mortality rates, but it disregards other sources of variance such as spatial variability. Thus: Then, the variance of the fraction of YOY present at the beginning of day d remaining at the end of the day (e −m d ) can be derived using the Delta method (Casella and Berger 2002): From this, we can calculate the variance of YOY t following an iterative process. The fraction of recruited YOY that remain alive at the end of day 2 (YOY 2 ) corresponds to: where m 1 and m 2 are the mortality rates on day 1 and 2 respectively. The variance of this product can be estimated using the usual approximate formula for the variance of two dependent variables (Goodman 1960): andρ is the estimated intra-class correlation for the cumulative remaining fraction of YOY. We have chosen to use the usual approximate formula, and not the exact formula, for the variance of two dependent variables because the additional terms included in the exact formula not present in the approximate formula incur into complex derivatives that, while representing substantial effort, add relatively little to the magnitude of the calculation (see Goodman 1960).
Similarly, the fraction of recruited YOY that remain alive at the end of day 3 (YOY 3 ) corresponds to: Fig. 1 Steps involved in the calculation of biomass at recruitment for fin-and shellfish species where derived life history tables exist, and productivity for species where derived life history tables are not available and its variance: Using this approach iteratively to YOY 4 , YOY t − 1 and finally to YOY t : Subsequently, the variance of D R can be calculated by applying to Eq. 13 the expression reported by Goodman (1960) that provides an unbiased estimate of the exact formula of the variance of the product of two independent random variables: where using the Delta method (Casella and Berger 2002): Finally, from the estimates of density at recruitment and their variances, we estimated a grand average of density at recruitment per habitat (D RH ) and its variance as: where D R is each of the estimates of density at recruitment for the given habitat, and N is the count of estimates in the habitat. Density estimates can be converted to biomass from knowledge of the mean individual fish weight at recruitment, which can be obtained from life history tables.
As indicated above, one of the main applications of our protocol is to allow for robust comparisons, or at least as robust as permitted by the size of the data sets available, of fin-and shellfish biomass across diverse shallow coastal habitats (i.e., SAV, oyster reefs, marshes, and non-vegetated bottoms). Such comparisons should be done for the same sampling time in all habitats; otherwise, temporal differences would confound the comparison and attribution of differences to habitat variability. If, rather than back-calculating to the time of recruitment, we had done comparisons across habitats with the same sampling time using the density data (or biomass after conversion from life history tables) directly reported in the data set, we would have been able to carry out only one comparison encompassing all habitats in the case of pinfish (i.e., month of May). All other comparisons for this species would have included three or four habitats, with different combinations of habitats for comparisons with the same number of habitats (for instance the comparison for July would encompass near non-vegetated, SAV, oyster reefs and marshes, and the comparison for October would encompass near non-vegetated, far non-vegetated, SAV and marshes, see Table 1). In the case of black drum, comparisons across habitats for the same sampling time would have involved at most three habitats, with many of them involving only two habitats (see Table 2).
The problems of reducing the number of habitat types that can be compared with the same sampling time when using directly reported density data, and additionally having discrepant combinations of habitat types among comparisons involving the same number of habitat types, apply to most other species included in the Hollweg et al. (2019) density data set. To avert these problems, we have developed the protocol presented above. The protocol allows for the simultaneous inclusion of all sampling times into the cross-habitat comparison by providing back-calculations from any sampling time to a common time point, i.e., the time of recruitment. The protocol provides an integrated and coherent comparison of fin-and shellfish biomass across habitats by bringing together all sampling times in the data set to the same time point. We have chosen time at recruitment because of its ecological and management significance (i.e., appearance of new recruits and onset for their growth in shallow coastal systems). Importantly, we propagate the error involved in our calculations throughout the derivation process, such that the final estimates allow for sound comparisons across habitats where the certainty and robustness of the differences found can be well informed.

Calculations of Productivity
Our protocol for the derivation of productivity values relies on estimation of the P:B ratio (ratio of productivity to biomass) and subsequent multiplication by biomass. Derivation of productivity using the P:B ratio and biomass has been carried out in the literature for macro-invertebrates (Sprung 1993, Cusson andBourget 2005) and fish (Waters 1977, Randall 2002). Here, we used an empirical model developed by Robertson (1979) for benthic macroinvertebrates that relates the species average P:B to its life span. The model corresponds to a linear regression fit using least-squares to the relationship between the base 10 logs of the two variables for 45 species of benthic macro-invertebrates including polychaetes, gastropods, bivalves, crustaceans, and echinoderms: In this equation, life span is expressed in years and ranges from 1.6 to 25.1 years. The P:B ratio is expressed in year −1 and ranges from 0.5 to 6.3 year −1 . It is important to stress that P:B corresponds to the mean ratio for the species, that is a ratio that includes all individual age classes and represents the mean P:B that would be measured at a population level including simultaneously all the different species' life cycle stages as they naturally occur in the system.
We can use the Robertson (1979) model to predict P:B from lifespan values for species of interest. Those predictions can be multiplied by mean biomass to derive estimates of productivity ( Fig. 1). Thus, the uncertainty in these productivity estimates comes from the uncertainty in the predicted P:B values and the uncertainty in the mean biomass. To derive the predicted P:B value, we inputted the lifespan for the species of interest in the equation above and solved for P:B. The entry corresponds to a species not included in the initial regression in Robertson (1979), and thus, the predicted P:B value for the entered lifespan is determined by the overall functional association between life span, body size and turnover rate across species (Brown et al. 2004), on the one hand, and idiosyncratic, species-specific variability in such functional association on the other. Therefore, we estimated the uncertainty of this prediction as the variance associated with a predicted single value of the dependent variable from a given value of the independent variable in the linear regression model (and not as the predicted mean value of the dependent variable for a given value of the independent variable, Neter et al. 1996).
This variance c varŶ single corresponds to: where MSE is the mean squared error from the model fit, n is the number of paired observations in the regression, x * is the specific value of the independent variable for which we seek the predicted value of the dependent variable, x is the mean for all the values of the independent variable used to obtain the regression fit, and S xx is the sum of squares of the independent variable. Since we used a regression model with base 10 log variables, the variance derived in this way corresponded to c varlog 10 P : B ð Þ: We used the Delta method (Casella and Berger 2002) to calculate c var P : B: where ln10 denotes the natural logarithm of 10, and P:B is the back transformed value of the predicted log 10 (P : B): To derive mean biomass and its variance by habitat, we first estimated the mean density and its variance in the specific habitat (Fig. 1). Estimates of mean density came from a two-step process. First, we pooled all density values corrected for overall gear efficiency corresponding to the same combination of habitat and sampling time, i.e., we derived D ht from D G ht , and c var D ht ð Þ from c var D G ht À Á as explained in Eqs. 9 and 10. Second, we estimated the grand average density (D H ) and its variance c var D H ð Þ ð Þ per habitat by pooling all the sampling times within the specific habitat: where N t is the number of sampling times in the habitat.
The purpose was to obtain a mean biomass value for the population that includes all individual age classes, or at least as many as possible, as they naturally occur. Thus, by first averaging all density values corrected for overall gear efficiency corresponding to the same combination of habitat and sampling time, this approach helps to reduce overweighting our final estimates with specific sampling times that could under-or overrepresent certain age classes. We then derived estimates of mean biomass per habitat (B H ) as the product between D H and mean individual weight (IW), with the latter also encompassing all size classes. Values of IWand its variance were obtained from the literature. We calculated the variance of B H using the equation reported by Goodman (1960) that provides an unbiased estimate of the exact formula of the variance of the product of two independent random variables: where values and variance of D H and B H are specific to the habitat, and the value and variance of IW are applied uniformly to all habitats. Finally, we derived estimates of mean productivity per habitat (P H ) as the product between B H and P:B, and its variance as (Goodman 1960): where the value and variance of the predicted P:B is applied uniformly to all habitats.

Species with Derived Life History Tables: Estimation of Biomass at Recruitment
Pinfish The data set used for pinfish in this paper comprised 278 density records from 26 papers (Appendix 3). Records per paper ranged from 2 (Reese et al. 2008, Cebrian et al. 2009) to 48 (Zimmerman et al. 1990). Density data spanned six habitats (marsh, oyster reefs, SAV, near non-vegetated bottom, far non-vegetated bottom, and near and far non-vegetated bottom combined together). Pinfish density data also spanned all twelve months and eight different gear types. Since the pinfish density data set compiled in Hollweg et al. (2019) was well populated with entries reported by month, we decided not to include entries reported by season (i.e., mean density values for a single or several seasons) in this paper to avoid the uncertainty associated with those entries when backcalculating density at recruitment. In total, the pinfish data set used here had 100 unique combinations of habitat, month, and gear type, with 51 of those combinations having one entry (Appendix 3). The number of entries per habitat/month/gear type combination ranged from 1 to 32 (marsh/May/drop sampler).
The density data compiled in Hollweg et al. (2019) for pinfish corresponds to YOY. Corrections for gear selectivity were not necessary since most YOY pinfish do not grow beyond 10 cm in length before they leave the shallow parts of coastal ecosystems to which they recruit (Nelson 1998) and, thus, generally remain below the maximum size threshold captured by the gear types included in the data set Minello 1997, 1998;see Hollweg et al. 2019 for further detail). Pinfish usually recruit to shallow areas through a few pulses over the winter and most YOY leave these areas towards deeper waters the next fall (Hansen 1969, Nelson 1998, Cebrian et al. 2009. Thus, to ensure that our calculations represent most of the recruited YOY, we assigned February 15th as the recruitment time for the back-calculations of density at recruitment from density values obtained at a later month, which we equated to the 15th day of the given month. In other words, day 1 was February 15th and day t was the 15th of the month at sampling for the estimation of density at recruitment from density values at a later date (Eqs. 12 and 13). Muncy (1984) reports that pinfish are about 20 days old when they recruit to shallow areas. Thus, the mortality rate assigned for February 15th was the mortality rate for 20-day-old fish, for February 16th it was the mortality rate for 21-day-old fish, and so forth consecutively through day 15th of the sampling month (Appendix 2). We only did these calculations for February through October, and not for November, December, and January since most YOY leave the shallow areas in the fall.
Pinfish biomass at recruitment (±SE) ranged from 0.007 ± 0.006 g DW per square meter in near non-vegetated habitat to 0.185 ± 0.060 g DW per square meter in oyster beds (Table 1, Fig. 2). In general, pinfish biomass at recruitment showed higher values in structured habitats, such as SAV beds, oyster beds, and marshes, than on bare sediment habitats, although recruitment values varied considerably within habitat types.

Black Drum
Density data availability for black drum in shallow coastal systems contrasted sharply with that for pinfish. Despite our extensive search (Hollweg et al. 2019), we could only find 37 records of black drum density from 4 papers (Appendix 3). The number of records per paper ranged from 3 (Gordon 2010) to 22 (Zimmerman and Minello 1984). The data set encompassed four habitat types (marsh, SAV, near nonvegetated bottom, and far non-vegetated bottom) and three gear types (drop sampler, throw trap, and seine). Due to the limited density data available for this species, we decided to use all the density data we could obtain and, unlike pinfish, we considered entries reported by month or season. Entries reported by season entail additional uncertainty in the backcalculation of density at recruitment. In total, there were 29 unique combinations of habitat, sampling time (month or season), and gear type, with only seven of them having more than one entry (Appendix 3). The number of entries per habitat/ sampling time/gear type combination ranged from 1 to 3 (near non-vegetated bottom, June, throw trap).
The density data compiled in Hollweg et al. (2019) for black drum corresponds to YOY. Corrections of this data for gear selectivity were necessary for enclosure-type gears (drop sampler and throw trap). This is because YOY black drum become larger than 10 cm before they leave the shallow areas in the fall (Sutter et al. 1986, Leard et al. 1993) and, thus, generally surpass the maximum size threshold for enclosuretype gears Minello 1997, 1998). Correction factors for selectivity can be derived from monthly individual size histograms of YOY fish. In this particular case, we can derive correction factors for specific months or seasons as the fraction of fish < 10 cm in the YOY population in the specific month or season. Despite our efforts, we could only find one report with such histograms for YOY black drum (Peters and McMichael 1990). The histograms are recreated in Fig. 3, and they show how YOY black drum grow quickly in shallow areas to reach a size over 10 cm by September. We calculated monthly correction factors as the average value of the two months in the two consecutive years in the study (1982 and 1983), and seasonal correction factors as the average for all months in the season pooling both years together (see legend of Fig. 3 for exact values). Since we derive these corrections factors based on a single report, we have chosen to not assign any variability to these factors (e.g., SE of the mean) but rather consider them as fixed constants for our calculations. The main purpose of this gear selectivity correction exercise is illustrative. Accurate derivation of selectivity correction factors will necessitate multiple reports of monthly individual size histograms.
Black drum YOY start recruiting to shallow areas in late winter and most recruitment has normally occurred by mid spring. The YOY stay in the shallow areas through the fall, at which time they move to deeper waters (Sutter et al. 1986, Peters and McMichael 1990, Leard et al. 1993). Based on this, for our calculations, we assigned April 15th as the recruitment time and included the months of April, May, June, July, August, September, and October, and the spring, summer, and fall seasons (with exceptions for enclosure-type gears, see Fig. 3). We excluded the months of November, December, January, February, and March, and the winter season. Day 1 was set at April 15th and day t corresponded to the 15th day of the month of sampling, or to the central point of the season of sampling (April 15th for spring; July 15th for summer; and October 15th for fall) for the estimation of density at recruitment from density values at a later date (Eqs. 12 and 13). According to Sutter et al. (1986), YOY black drum are around 15 days old when they arrive in shallow areas. Thus, the mortality rate at day 1 corresponded to the mortality rate for 15-day-old fish, at day 2, it corresponded to the mortality rate for 16-day-old fish, and so forth consecutively through day t (Appendix 2).
Biomass at recruitment was much lower for black drum than for pinfish in the habitats studied (Table 2, Fig. 2). Values (± SE) ranged from a rounded-up value of 0 by the authors (± 0.015) in far non-vegetated, to 2.4 × 10 −5 ± 0.009 in SAV beds, to 0.003 ± 0.005 in marshes, and to 0.011 ± 0.008 g DW per square meter in near non-vegetated habitat.

Mud Crab
The density data set compiled by Hollweg et al. (2019) for this species contained 162 records from 19 papers (Appendix 3). The number of records per paper ranged from 1 (Zeug et al. 2007, Roth andBaltz 2009) to 46 (Zimmerman et al. 1990). Density data spanned all six habitats and three gear types. Data for mud crab were reported by single month, by combinations of two months within a season, by single seasons, or even by combination of seasons. We included all these sampling times in our calculations since the intent for this species is to generate density estimates that integrate all individual age classes as they naturally occur in the habitats studied. In total, there were 49 unique combinations of habitat, sampling time, and gear type, with 17 of those combinations having one entry (Appendix 3). The number of entries per habitat/sampling time/gear type ranged from 1 to 20 (marsh/May/drop sampler).
Corrections for gear selectivity were not necessary since mud crabs do not normally exceed 10 cm in length (Williams 1984). Most of the sampling times averaged within each habitat to calculate the grand average density in the habitat corresponded to spring, summer, fall, or combinations of those, with little sampling done in winter. Hence, since winter is poorly reported, our final estimates of grand average density per habitat may have some bias in relation to well-balanced estimates that encompass all of the species age classes as they occur in nature. However, data were reported in a similar frequency for spring, summer and fall within each of the habitats. In other words, for each habitat, approximately 1/3 of the data corresponded to spring, 1/3 to summer, and 1/3 to fall. Therefore, despite potentially having some implicit bias, our estimates should still provide sound comparisons across habitats because spring, summer and fall are sampled to a similar extent (1/3 of the total data for the habitat) within each of the habitats.
Mud crab productivity (±SE) in the habitats studied ranged from 0.180 ± 0.569 (far non-vegetated) to 0.793 ± 2.516 (near non-vegetated) in non-structured habitats, and from 0.425 ± 1.388 (marsh) to 4.098 ± 5.567 g DW per square meter per year (SAV beds) in structured habitats (Table 3, Fig. 2). In general, productivity values were higher in structured than in non-structured habitats. However, large variability in mud crab productivity was observed within each of the habitats studied.

Gulf Stone Crab
Data availability was much more limited for Gulf stone crab in relation to mud crab. Despite our extensive search (Hollweg et al. 2019), we only obtained 34 density records from 6 papers for Gulf stone crab (Appendix 3). The number of records per paper ranged from 1 (Roth and Baltz 2009) to 18 (Peterson and Stricklin 2008). Density data spanned five habitats and four gear types. Data for Gulf stone crab was reported by single months, combinations of months, seasons, and combinations of seasons. Similarly to mud crab, we included all these sampling times in our calculations since we intended to generate density estimates that encompass all individual age classes as they naturally occur in the habitats studied. In total, there were 20 unique combinations of habitat, sampling time, and gear type, with 13 of those combinations having one  Peters and McMichael 1990). The dashed green light shows the fraction of fish caught with enclosure-type gears (drop sampler and throw trap) assuming a maximum size threshold of 10 cm Minello 1997, 1998). These fractions have been applied for gear selectivity corrections, and correspond to 0.65 for July, 0.165 for August, and 0.605 for summer (see text for more details). Records for September, October, and the fall obtained with enclosure-type gears are excluded in our calculations since very few YOY would be captured with these gears at those sampling times according to the histograms Table 3 Mud crab: Density at sampling and SE; overall density at sampling per habitat and SE; and productivity and SE. Density is in number of individuals per square meter of habitat, and productivity is in gram DW per square meter of habitat per year. Productivity has been derived by multiplying overall density times the mean individual dry weight ± SE (0.563 ± 0.014 g per individual) and production to biomass ratio ± SE (1.841 ± 0.128 per year) encompassing all age classes. In turn, the mean production to biomass ratio is calculated from the species life span of 3 years using Robertson (1979). Sources are Turoboyski (1973) for mean individual weight, and Williams (1984) Table 4 Gulf stone crab: Density at sampling and SE; overall density at sampling per habitat and SE; and productivity and SE. Density is in number of individuals per square meter of habitat, and productivity is in gram DW per square meter of habitat per year. Productivity has been derived by multiplying overall density times the mean individual dry weight ± SE (3.768 ± 0.636 g per individual) and production to biomass ratio ± SE (1.059 ± 0.074 per year) encompassing all age classes. In turn, the mean production to biomass ratio is calculated from the species life span of 7.5 years using Robertson (1979). Sources are Caudill (2005) and Ricciardi and Bourget (1998) for mean individual weight, and Gerhart and Bert (2008)  Since Gulf stone crabs do not normally exceed 10 cm (Gerhart and Bert 2008), corrections for gear selectivity were not necessary. Except for near non-vegetated bottom, sampling times were similarly distributed throughout the year within the habitats examined, with similar relative sampling frequencies reported for spring, summer, fall, and winter within each habitat. This suggests our density estimates for Gulf stone crab should encompass all individual age classes as they naturally occur in the habitats examined, except in near nonvegetated bottom.
Values (± SE) of Gulf stone crab productivity showed large variability within and across habitats, ranging from 0.004 ± 0.846 in near and far non-vegetated, 0.018 ± 1.117 in SAV beds, 0.617 ± 1.882 in marshes, 1.141 ± 2.238 in oyster beds, and 3.108 ± 4.510 g DW per square meter per year in near non-vegetated habitats (Table 4, Fig. 2).

Discussion
In this paper, we combine a number of established methods to develop procedures that generate estimates of biomass at recruitment for species of fin-and shellfish with derived life history tables, and estimates of productivity for species where sufficient information does not exist to derive life history tables. Our results contribute to a growing battery of tools for the calculation of fish biomass and productivity in coastal ecosystems (e.g., Wong et al. 2011, French McCay et al. 2003a,b, Nelson et al. 2013, French McCay et al. 2015. The calculations presented here are straightforward and can be implemented in guidelines and policies for environmental management. For instance, our protocol can help managers estimate the biomass and productivity of species targeted for commercial, conservation, or restoration purposes. It also provides the base to estimate fisheries productivity under various scenarios of conservation and/or restoration of coastal habitats such as SAV beds, marshes, and oyster reefs. Importantly, the protocol further quantifies the uncertainty around all such estimates, thereby helping managers make well-informed decisions. Our calculations indicate that pinfish biomass at recruitment is higher in structured (SAV beds, oyster reefs and marshes) than in non-structured habitats (bare sediment bottoms). These results are consistent with the preference that recruiting juvenile pinfish typically show for structured vs. non-structured habitats. Indeed, other studies have also shown higher levels of recruitment, abundance, and growth of juvenile pinfish in structured vs. non-structured habitats (Stoner 1983, Jordan et al. 1996, Tolan et al. 1997, Cebrian et al. 2009. Black drum biomass at recruitment is generally much lower than for pinfish in the habitats studied. In addition, our results suggest that recruiting juvenile black drum does not show the preference for structured habitat that is apparent with juvenile pinfish, which is consistent with past reports on the life cycle, biology, and habitat distribution for black drum (Osburn and Matlock 1984, Cody et al. 1985, Sutter et al. 1986). At any rate, density data for juvenile black drum are very scarce for the habitats studied; thus, the results must be regarded with caution. Because of the low sample size, it is likely that the estimates of biomass at recruitment for this species are inaccurate (i.e., rather far from the true value) and imprecise (i.e., high SE in relation to the mean value). Thus, while it seems possible that biomass at recruitment is generally lower for black drum than for pinfish in the habitats studied, our calculations cannot resolve with accuracy how much lower it is and also the true differences (or lack of) in black drum biomass at recruitment across the habitats studied. New measurements and resulting larger data sets of YOY density for this species are needed to analyze these differences with rigor.
Our estimates of mud crab productivity are generally larger in structured than in non-structured habitats. However, there exists large variability in these estimates within any given habitat. Indeed, the ratio of SE to the mean value of productivity in the habitat varies from 1.4 in SAV beds to 2.1 in oyster beds to around 3.5 for the rest of the habitats. This high variability emerges, at least in part, from the high variability also observed in the primary density data gathered for our calculations. For instance, when comparing the means versus their SE's of density values corrected for overall gear efficiency corresponding to the same combination of habitat, sampling time, and gear type (D G ht ), we find a lower value for the mean than for its SE in 46 out of the 49 combinations (higher mean than its SE in 3 of the combinations). In contrast for pinfish, we find a lower value for the mean than for its SE in 52 out of the 100 combinations (higher mean than its SE in 43 combinations, and equal mean and SE in 5 combinations). This higher variance in mud crab in relation to pinfish density is not related to the number of records averaged in specific combinations of habitat, sampling time, and gear type when there are multiple records for the combination, since overall those numbers were similar in the two species (see Appendix 3). The differences in density variability between mud crab and pinfish pervade as averages are calculated across gears for specific combinations of habitat and sampling time ( D ht ), although somewhat muffled. In addition, the high variability in density values for mud crab in relation to pinfish seems irrespective of the length covered by the sampling time, since mud crab density values show similarly large variability when comparing sampling times including one month, several months, and several seasons. Thus, when comparing mud crab and pinfish, the two relatively well documented species examined here, it appears that mud crab density is intrinsically highly variable in relation to pinfish density (see also Boyle et al. 2010, Gagnon and Boström 2016). This substantial variability in the compiled density values cascades through the final productivity estimates, and it hinders strong inferences regarding differences in mud crab productivity across the habitats studied. Larger data sets for this species in shallow coastal systems may alleviate these caveats by decreasing the variability around the estimates; however, our search has been exhaustive and only, if any, little additional data should exist at present.
Our results for Gulf stone crab show high variability in its productivity within and across the habitats studied. The calculations are based on a poorly reported data set, and thus, these results must be viewed with caution. As a matter of fact, the estimate for near non-vegetated habitat comes from only one mean density value for a sole sampling time, and all other estimates also suffer from highly limited sample size. Because of this, it is likely that our estimates of Gulf stone crab productivity are inaccurate (i.e., they probably are considerably off from the true value), and also imprecise (as seen in the magnitude of the SE in relation to the mean value). Thus, our calculations may offer a general idea of Gulf stone crab productivity in coastal habitats, but they bear little value for strong inferences regarding the accurate assessment of this productivity and its differences across the habitats studied. Because our compilation is exhaustive, it seems this is the best information our calculations can provide until more density data are obtained for this species in shallow coastal systems.
Despite our efforts to estimate the variability around the derived values of biomass at recruitment and productivity, there are sources of variability not included in our calculations. For instance, we do not include the variability across different populations corresponding to the same combination of habitat, sampling time, and gear type, nor do we include the variability due to author bias that can be generated when several density values for the same combination come from the same authors or paper. Additionally, our estimates of fraction of recruited YOY remaining at a given sampling time do not include the spatial variability of the daily YOY mortality rates derived in life history tables (we only include temporal variability using Bradford (1992)). We made an exhaustive search but, unfortunately, we did not find any reports that would have allowed us to derive robust estimates of such variability and include it in our calculations. The best information we found was measurements of YOY mortality rates obtained from the decrease in YOY abundance over a month's time in spring reported by Nelson (1998) in Choctawhatchee Bay, Tampa Bay, and Charlotte Harbor. Given the closeness of the values obtained (0.022, 0.021, and 0.023 day −1 , respectively) and locations studied, we did not judge these values as appropriate to derive robust estimates of spatial variability in YOY mortality rates.
Our calculations also include not well-constrained sources of uncertainty. Based on similarities between gear types and habitats, we assigned surrogate correction factors for gear type and habitat combinations where the literature lacked specific correction values, which may generate bias in our estimates of biomass at recruitment and productivity (see detailed explanation in Hollweg et al. 2019). The inclusion of seasons along with months as sampling times in the calculations of density at recruitment for black drum may generate further bias. Including seasons along with months increases the sample size; however, this also generates an additional source of uncertainty in that the sample collection date is assigned to the central point of the season (April 15th in spring; July 15th in summer, October 15th in fall), but the real collection could in reality have happened any time during the season and be substantially off from the assigned collection date. Ultimately, the choice of including seasons as sampling times in the calculations should be dictated by whether, despite the uncertainty brought about by assigning the collection date to the central seasonal point, the derived estimates are more robust than if seasons are not included. There are quantitative tools to help inform such decision, such as error propagation techniques and sensitivity analysis (Lehrter and Cebrian 2010). Here, we have not carried out analyses to justify inclusion of seasons along with months as sampling times in our calculations. We decided to include seasons based on the fact that about 40% of the mean density values corrected for overall gear efficiency used to back-calculate density at recruitment from density obtained at a later date (D ht ) were reported by season, so discarding those values would have further crippled an already highly limited data set. However, had we elected to exclude them, we would have reached similar results: black drum recruitment seems lower than for pinfish and not associated with structured habitats, but a larger data set is needed to substantiate these suggestions.
Therefore, large density data sets simultaneously available for YOY pinfish and black drum, and for total populations of mud crab and Gulf stone crab, along with more complete gear efficiency data sets in shallow coastal systems would much help produce accurate (close to the true value), precise (with relatively low SE in relation to the mean value), and inclusive (including as many sources of variability as possible) estimates of biomass at recruitment for pinfish and black drum and productivity for mud crab and Gulf stone crab. Such large data sets would allow for robust comparisons and inferences of differences among species and habitats. Our work also highlights the importance of documenting data sets thoroughly. When reporting information, an effort should be made to provide the most complete documentation possible regarding means, sample sizes, effect sizes, and variability (SD and SE), thereby allowing for rigorous further analysis of the data. For instance, when using our imputation method to estimate missing SE's, incomplete data reporting (e.g., no sample size reported) generates bias and higher variability for the density weighted average derived for specific combinations of habitat, sampling time, and gear type (Hollweg et al. 2019). Large data sets that are well documented can do much in advancing our understanding of fin and shellfish biomass and productivity in coastal ecosystems as well as our capacity to manage these ecosystems.
The results for pinfish and mud crab point to differences across the habitats examined in biomass at recruitment for the former species and productivity for the latter, suggesting higher values in structured vs. non-structured habitats. However, when attempting to assign statistical significance to these visually apparent differences, our efforts prove challenging. This is an important endeavor since the extent of within-habitat variability found for the two species sheds doubt on whether statistically significant differences do occur among habitats, particularly for mud crab productivity. In our calculations, the density values corrected for overall gear efficiency corresponding to the same combination of habitat, sampling time, and gear type (D G ht ) have a non-normal distribution. When pooling these values into the calculation of D ht and D H for mud crab, these summatory variables are based on too few summands to, based on the Central Limit Theorem, reasonably assume they have a normal distribution. Nonnormality pervades through the final calculations of B H and P H . For pinfish, the calculated values of density at recruitment (D R ) have a non-normal distribution, and similarly the summatory variable D RH derived from D R is based on too few summands to assume normality. Biomass at recruitment per habitat should also have a non-normal distribution, since it is derived as the product between D RH and a constant. Thus, we could not simply compare estimates of biomass at recruitment for pinfish, or productivity for mud crab, among habitats using a z-statistic. In the same vein, because we cannot assume normality for many of the variables derived throughout our calculations, it is highly questionable to use the Welch-Satterthwaite formula (Ku 1966, Lehrter andCebrian 2010) to meticulously propagate degrees of freedom with the intent of deriving confidence intervals for the final estimates of biomass at recruitment or productivity per habitat based on the t-statistic.
The possibility of using bootstrapping methods to build confidence intervals for final estimates of biomass at recruitment and productivity per habitat is also questionable within our calculation framework, since we end up with a rather limited number of estimates of density at recruitment for pinfish (D R , from 2 estimates in oyster reefs to 9 in various other habitats) or density at sampling for mud crab within habitats (D ht , from 4 in far non-vegetated and oyster reefs to 12 in near non-vegetated and marshes) and the confidence limits derived in this way would not be very robust. In the case of pinfish, perhaps bootstrapping techniques could be applied more powerfully if estimates of density at recruitment were derived from each primary single record, that is for every single record pooled into density weighted averages for specific combinations of habitat, sampling time, and gear type when several records exist for the specific combination, and not from pooled averages for the same sampling time (D ht ). This, however, necessitates further analysis since the approach would be significantly different from the one used here. Regarding mud crab, density values need to be averaged by sampling time to prevent over-representation of a given sampling time when deriving productivity estimates, so the use of bootstrapping seems inherently limited when using the P:B method. We suggest the use of Bayesian methods to test for inferences regarding differences in our estimates of biomass at recruitment and productivity among species and habitats. Given reasonably large density data sets, as well as sufficient detail for other complementary information such as thorough gear efficiency correction factors, Bayesian methods should prove adequate to derive sound conclusions regarding differences in biomass at recruitment and productivity between species and habitats as calculated with our approach, and we suggest this as a promising avenue of research for a more complete utilization and application of our approach.
Our approach can be extended to derive productivity estimates for fin-and shellfish species with derived life history tables. For instance, we can adapt our calculations to estimate productivity of pinfish and black drum from hatching to a later time in their juvenile life stage before migration to deeper waters. To do this, we would first calculate the density values corrected for overall gear efficiency for each of the sampling times considered (D ht ), be months or seasons, and from this we would calculate the density at the time of interest. If the time of interest is earlier than the sampling time for D ht , we would then back-calculate density at the time of interest from D ht . These back-calculations would be tantamount to the calculations done here (Eq. 13), only that the time of interest would not be time at recruitment. If the time of interest is later than the sampling time for D ht , we would fore-calculate density at the time of interest from D ht . Fore-calculations would follow the same rationale as back-calculations, only that we would be moving forward in time, rather than backwards. In this regard, Eq. 13 can be modified to calculate density at a later time from density at a former time. Once we have all estimates of density at the time of interest (i.e., one independent estimate of density at the time of interest from each D ht ), we would calculate their grand average and multiply times the mean individual fish weight at that time, which can be obtained from the derived life history tables. These final estimates would correspond to fish productivity from hatching to the time of interest, and it would be expressed in g DW per square meter per the time period elapsed from hatching to the time of interest. The variability (SE) for these productivity estimates would also be obtained following the procedure presented here, conveniently adapted to the particular amendments made. Similarly, we could estimate juvenile biomass export by equating the time of migration to deeper waters to the time of interest in the calculations. As a matter of fact, our approach could be used to estimate fish productivity from hatching to any life-cycle stage (age class), assuming that there are good age-class-specific density data for the species.

Conclusion
Building from existing methods, here we developed a protocol that allows for calculations of biomass at recruitment for species of fin-and shellfish with derived life history tables, and for calculations of productivity for species where derived life history tables do not exist. Our procedure propagates the uncertainty of variables included in the calculations, and thus generates brackets of variability for the final estimates. For species with better reported density data sets, reasonable inferences regarding differences among species or habitats can be suggested. For species with poorly reported density data, only a general idea of biomass at recruitment or productivity can be gleaned. However, larger data sets of fin-and shellfish density in shallow coastal systems, and other complementary information such as thorough gear efficiency correction factors, are needed to generate estimates of biomass at recruitment and productivity that are more accurate, precise, and inclusive of variability sources. In combination with larger data sets and Bayesian statistics, our protocol offers promise for a better understanding of fisheries productivity in coastal systems and enhanced management of these systems.
Acknowledgments We thank Josh Goff for the help with table and figure preparation. The ideas presented here benefited from discussions with Charles Peterson, Sean Powers, and Lawrence Rozas, whom we also wish to acknowledge. Jennifer Doerr, of the NOAA Southeast Fisheries Science Center, and Tony Marshak, of ECS Federal LLC in support of the NOAA Fisheries Office of Science and Technology, reviewed the manuscript and provided useful comments. We also want to acknowledge two anonymous reviewers and the editors of Estuaries and Coast for their useful input during the review process. The scientific results and conclusion of this publication, as well as any views or opinions expressed herein, do not necessarily represent the view of the other natural resource Trustees for the BP/Deepwater Horizon NRDA. Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government.

Appendix 1
Equations used for the derivation of individual fish mass and mortality rates. The equations are run on daily steps for the first year of life (see Appendix 2). Species-specific values for initialization and validation of the equations have been obtained from Nelson (2002) for pinfish and Murphy and Muller (1995) for black drum. Water temperature was set at 26°C, which is a representative mean value during the larval stage period for these two species in the region studied (Boyer et al. 2011). For details, see pages 11-13 in French McCay et al. (2015) Larval stage: Duration of the larval stage (Houde 1989) where D l is the duration of the larval stage (days) and T is water temperature (C) Growth/Length (Pepin 1991) where G l corresponds to the daily growth increment of the larval individual (mm/day), T is water temperature (C), and L is larval length (mm) Mortality (Pepin 1991) where M l is the larval mortality rate (day -1 ), T is water temperature (C), and L is larval length ( where G j is the juvenile growth rate (day -1 ) from the end of the larval stage to age 1 year, W 1 is the individual wet weight at age 1 year (grams), W l is the individual wet weight at the end of the larval stage (grams), D j is the duration of the juvenile stage up to age 1 year (days), W τ is the individual wet weight at age τ (grams), and D l is the duration of the larval stage (days).
Mortality (first year of life) (Lorenzen 1996) where M j is the juvenile mortality rate (day -1 ), and W is the individual wet weight (grams) Length-weight conversion equations Length to wet weight (Wiebe and Davis 1985) Wet weight to dry weight (Nixon and Oviatt 1973) where W wet is wet weight (grams), L is length (mm), and W dry is dry weight (grams)

Appendix 2
Estimates of mortality rates and individual fresh weight for daily intervals during the first year of life for pinfish larvae and juveniles. For calculation details, see Appendix 1. Estimates of mortality rates and individual fresh weight for daily intervals during the first year of life for black drum larvae and juveniles. For calculation details, see Appendix 1 (Tables 5 and 6)                  Ballestero, M.T. Huisenga, and K.G. Benson. 2019. Metaanalysis of nekton utilization of coastal habitats in the northern Gulf of Mexico. Estuaries and Coasts (in press). Unique combinations of habitat, sampling time, and gear type for the four species studied here. The number of records for the combinations is also provided (Tables 7 and 8)    Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.