Frequency, Magnitude, and Possible Causes of Stranding and Mass-Mortality Events of the Beach Clam Tivela mactroides (Bivalvia: Veneridae)

Stranding combined with mass-mortality events of sandy-beach organisms is a frequent but little-understood phenomenon, which is generally studied based on discrete episodes. The frequency, magnitude, and possible causes of stranding and mass-mortality events of the trigonal clam Tivela mactroides were assessed based on censuses of stranded individuals, every four days from September 2007 through December 2008, in Caraguatatuba Bay, southeastern Brazil. Stranded clams were classified as dying (closed valves did not open when forced) or dead (closed valves were easily opened). Information on wave parameters and the living intertidal clam population was used to assess possible causes of stranding. This fine-scale monitoring showed that stranding occurred widely along the shore and year-round, with peaks interspersed with periods of low or no mortality. Dead clams showed higher mean density than dying individuals, but a lower mean shell length, attributed to a higher tolerance to desiccation of larger individuals. Wave height had a significant negative relationship to the density of dying individuals, presumed to be due to the accretive nature of low-energy waves: when digging out, clams would be more prone to be carried upward and unable to return; while larger waves, breaking farther from the beach and with a stronger backwash, would prevent stranding in the uppermost areas. This ecological finding highlights the need for refined temporal studies on mortality events, in order to understand them more clearly. Last, the similar size structure of stranded clams and the living population indicated that the stranded individuals are from the intertidal or shallow subtidal zone, and reinforces the ecological and behavioral components of this process, which have important ecological and socioeconomic implications for the management of this population.


Introduction
Stranding and mass mortality on the shore are common events for many invertebrate groups, and their frequency, extent, and duration have been increasingly reported in recent decades [1][2][3]. These events are attributed to several factors, such as storms, positive sea-temperature anomalies, parasitism, harmful algal blooms, and hypoxia [1,2,[4][5][6][7][8][9][10][11][12]. Despite the ecological importance of these processes, data on the frequency and magnitude of stranding and massmortality events of invertebrates are generally limited in time, probably because of the sudden and generally episodic occurrence of this phenomenon [2,3,13], for sandy-beach clams in particular [2,9,10,[13][14][15][16][17]. Mass mortality of sandy-beach clams has important ecological implications because it can affect the density of these frequently dominant populations, altering the structure of the benthic community [18,19], as well as the organisms that rely on this food source [20]. In extreme cases, mass mortality can decimate clam populations along their geographical range, as reported for the yellow clam Mesodesma mactroides and the surf clam M. donacium [9,11,[21][22][23]. In addition, in view of the high importance of sandy-beach clams in recreational and artisanal fisheries of developing countries in South America, mass-mortality events of these organisms can be socioeconomically significant [18,22,[24][25][26]. The trigonal clam Tivela mactroides is an important sandy-beach species, both economically and ecologically [20,[27][28][29][30][31][32][33][34]. It has a wide distribution from Venezuela to Brazil (Pará to Santa Catarina states) and in the West Indies [35]. This suspension-feeder inhabits the intertidal and subtidal zones of sandy beaches, where it is often the dominant macrofaunal component [28,[36][37][38]. The trigonal clam frequently reaches high abundance in Venezuelan beaches, where it is an important resource of artisanal fisheries [28,30,31]. On these tropical beaches, no studies of mass mortality of the trigonal clam have been conducted, although there are brief notes on observations made under conditions of strong winds and waves [28]. The population of the trigonal clam living in the intertidal zone of sandy beaches in Caraguatatuba Bay, southeastern Brazil, has continuous recruitment, high individual growth and mortality rates, and reaches high levels of density, biomass and secondary production [32,33]. On several occasions, large numbers of stranded clams were observed in the intertidal zone of these beaches, and these mass strandings were presumed to be related to population density-dependent processes such as individual growth and mortality rates [32]. Also for this population, an assessment of harvesting intensity showed that under ordinary conditions, the amount of clams harvested yearly is about 25 tons. The same study showed that the mass-mortality events are also socioeconomically important, since a decrease in harvesting was attributed partly to these events [34]. The authors suggested that the foul odor and unpleasant sight caused by the stranded and dead clams were responsible for a decrease in recreational use of the beaches in Caraguatatuba Bay. Mass-mortality events can have dramatic effects on fishing and management of natural resources [11,22]. Moreover, some physical changes related to global climate change can increase the frequency and severity of mass-mortality episodes on sandy beaches [9,11]. In view of these considerations, assessment of mortality events of T. mactroides is a crucial step in managing populations of this species.
Among the many possible causes of mortality are aspects of the population dynamics such as overpopulation, competition for resources, and parasitism, but also physical features, since these clams live under harsh conditions, depending largely on the characteristics of the sediment where they burrow and on the wave energy across the shore [9,11,17,32]. Wave energy has a high potential to transport sediment as well as individual clams, while the sediment conditions may favor or preclude their burrowing and survival [39]. In fact, shoreward passive migration of small-sized clams was hypothesized to explain spatial partitioning between small (subtidal) and large (intertidal and shallow subtidal) individuals [37].
The present study aimed to shed some light on stranding and mass-mortality events of these beach clams, by assessing their frequency, magnitude and possible causes for the population of T. mactroides in Caraguatatuba Bay. A comprehensive study was carried out on stranded individuals, both dying and dead, at a 4-day temporal scale during 15 months. The concomitant physical and biological conditions were assessed in relation to the stranding events. Given the general importance of hydrodynamics for mortality events on sandy beaches [40], wave parameters were considered in order to test the hypothesis that these are related to the number or size of stranded individuals. Similarly, the data for intertidal-living buried individuals, obtained by Turra et al. [32], were also included as a potential factor influencing the structure of stranded individuals, i.e., if the density or size of the stranded clams is related to the respective parameters of the buried intertidal individuals.

Study area and sampling
The coastline of Caraguatatuba Bay (23°37'S to 23°44'S and 45°24'W to 45°26'W) in southeastern Brazil is about 16 km long and is uninterrupted by physical barriers. The sandy beaches in the study area have a gentle slope (1/33-1/41 m), a wide surf zone (ca 150 m), fine to very fine sand (0.8-1.12 mm), and waves predominantly with a period of 6 s and height of 0.4 m [32,33,37].
The mortality of Tivela mactroides was assessed in a 4-km area of the central-southern portion of the bay. This area was selected because intertidal densities of this clam are higher in this part of the bay [32,37]. Sampling was carried out under license granted by the Instituto Brasileiro do Meio Ambiente e dos Recursos Naturais Renováveis (IBAMA-DIREN No. 08/2001), at four-day intervals from September 8, 2007 to December 3, 2008. Preliminary evaluations indicated that this periodicity was sufficient to include all mortality events in the samples. Individuals stranded on the beaches were sampled in three subareas, each 400 m long (parallel to the waterline) and 15 m wide (perpendicular to the waterline), which were randomly selected on each sampling day in the study area (Fig 1). This 15-m strip was established 30 m below the vegetation edge (supralittoral limit), which was not included in the sampling, in order to avoid the effect of debris removal by the local government. Also, the area below the 15-m strip, which corresponded to the lower intertidal zone, was not included in the sampling design due to the impossibility of sampling during neap low tides. In contrast to the higher levels, where the stranded individuals generally had open valves, in this 15-m strip most of the stranded clams were found with the valves closed, indicating that they reached the beach strip recently.
In each subarea of 400 x 15 m, a rectangular grid was then defined, in which a sampling station was systematically established every 20 m (0, 20, 40, 60, . . . 400 m) along the longitudinal axis of the subarea, totaling 21 sampling stations (Fig 1). Across the shore, a number among 15 possibilities (0, 1, 2, 3,. . .15 m) was randomly selected for each sampling station. The 21 sampling units in each subarea were examined using a 1 m 2 quadrat. In each quadrat, all stranded individuals of T. mactroides were sampled and separated into dying (closed valves that did not open when forced) or dead (open valves or easily opened when forced). The clams were placed in labeled plastic bags, and subsequently counted and measured for total length (in mm).

Wave modeling
The height, period and power of the local waves were estimated daily for the entire sampling period (September 1, 2007 to December 3, 2008). Based on the regional offshore wave information extracted from the global wave generation model WaveWatch III (NCEP/NOAA) [41], a numerical model was applied in order to propagate the waves onshore until they reached the area of interest. The numerical model was the Delft3D (open-source version from Deltares) WAVE module, which simulates the propagation and transformation of wind-generated waves as they move over varying bottom morphology [42]. It simulates effects such as wave refraction, diffraction, shoaling and set-up in coastal areas.
The Caraguatatuba model used in this application covered most of the northern coast and continental shelf of São Paulo State, keeping the area of interest well inside the model domain. Grid resolution varied in the domain, with increased resolution onshore. Boundary wave data, defined at the shelf break, represented the wave climate (significant wave height and corresponding wave period) for the entire period of analysis, with information given at six-hour intervals. Subsequently, near-shore wave data were obtained at the entrance of Caraguatatuba Bay (23°40'S and 45°22'W), providing the local wave conditions during the sampling period. The daily mean wave height and period were calculated from four daily values. In order to account for the synergy of wave height and period, the wave power P was estimated by: P = ρ g 2HT/32Π, where ρ is water density (1,027 kg/m 3 ), g the acceleration due to gravity (9.81 m/s 2 ), H the wave height (m), and T the wave period (s). P is given in W/m.

Data on living buried individuals
Considering that some population parameters (e.g., density and shell length) of the intertidal buried individuals can provide important information on the stranding process, the data for these clams obtained by Turra et al. [32] were also used. This study sampled buried individuals in summer and winter 2007 and 2008. Each sampling period included three random intervals of spring low tides, when living, buried clams were collected, counted and measured for total length (in mm). Data analysis General characterization. The mean density (ind.m -2 ; total number of occurrences divided by 21, which is the mean number of individuals per 1 m 2 quadrat within the area) and mean shell size (mm) were calculated separately for dying and dead stranded individuals for each area sampled. The strength of the relationship between dying and dead clams was evaluated for both density and shell size (for all cases the ordinary least square regression was appropriate).
Another goal was to assess temporal differences, as well as if there were differences according to the individuals' condition (dying or dead), in mean density and size of stranded individuals. Therefore, for each biotic measurement (density and size), months and condition were considered as factors in a factorial ANOVA for mixed design, where data on dead and dying individuals within a sample were considered repeated measures. The assumptions were all assessed and met, except for density, which showed a wide range of values, and consequently a high heterogeneity of the variances. Therefore, the density data were transformed (log 10 (x+1)) to meet the ANOVA requirements [43].
Relationship to wave conditions and to the living population structure. To assess whether the wave regime was related to the stranding of T. mactroides, mean daily values for density and size of stranded individuals were modeled separately according to wave height, period and power. Because the population may show a delayed response to a change in environmental conditions, for the wave data we used the parameters from the sampling day (time 0) and from the preceding days, using a time lag from one (time -1) to four days (time -4). Because these wave data are highly correlated, exploratory analysis and nested model selection were applied to determine which parameters should be used in the model without the inclusion of correlated factors. Model selection was performed [44,45] using a quantile regression analysis, more precisely the upper limit [46,47]. This technique shows the extreme state of a relationship, here considered of key importance because other factors, for example the lack of living individuals, would result in values well below this limit, and not contribute to the true shape of the relationship between the variables [48].
This modeling process indicated that the use of dying clams alone responded better to the environmental factors, and therefore better than the total number of stranded clams. Therefore, data on dead individuals were not used in model analyses. Because the transects were placed randomly in each sampling event, and the total area was too large to be fully cleared between sampling events, dead clams could have been pseudoreplicated. However, although among the stranded individuals a dead clam could be present longer than 4 days, a dying clam had certainly stranded recently enough to prevent an overlap between sampling periods.
Optimum model selection was performed by applying backward stepwise selection based on Akaike criteria and using hypotheses tests for each model run. Graphical analysis of the residuals was used to check for homogeneity of variances, normality and independence, until an ideal model meeting all assumptions was met [44]. Here, because a clear non-linear relationship was found, a generalized additive model (GAM), with Gaussian distribution and using a regression spline (thin plate) as a smoothing basis, was applied using the mgcv R package [49].
Data from Turra et al. [32] on the structure of the intertidal stratum of the living population were also used to investigate a possible relationship of this part of the living population to the structure of the stranded individuals. These data were sparser in time than the present study, and therefore only the corresponding monthly values were isolated and used in the analysis, which included all stranded individuals. Similarly to the previous analysis, exploratory data assessments were performed to validate this choice, and as expected given the broader time scale, the use of both dying and dead individuals proved to be better. The shell length and density of living individuals, obtained from Turra et al. [32], were compared with the shell length and density of the present stranded individuals, and an ordinary least square regression was used to assess the relationship in both cases, followed by a paired t-test for shell length.

General characterization
A total of 118,871 stranded clams, of which 23,424 were dying (20%) and 95,447 were dead (80%), were collected in the 75 sampling events. Most of the quadrats and subareas sampled contained both dying and dead stranded individuals. Relationships between dying and dead clams were positive for both density (F 1,223 = 253.5, p<0.01, r 2 = 0.53) and shell length (F 1,166 = 105.6, p<0.01, r 2 = 0.39), showing that there was a significant, though weak, correspondence between dead and dying clams. The density of stranded clams differed significantly over time (F 15,32 = 4.40, p<0.01), with an overall mean value of 27.09 ind.m -2 (± 45.87 SD) and wide oscillations throughout the study period (Fig 2A). Clam condition also differed significantly over time ( In 2008, high densities of stranded individuals were observed from January through June, with oscillating frequencies, and then decreased, with isolated peaks between August and October. In the last sampling events, stranded individuals were recorded only once. Regarding mean shell length, dead clams were recurrently smaller then dying ones (F 1,26 = 87.70, p<0.01; 17.06 ± 2.16 mm and 18.37 ± 2.46 mm, respectively) ( Fig 2B). There was also a significant temporal difference in shell length (F 12,28 = 37.41, p<0.01), with no interaction with individual condition (F 1,26 = 0.75, p = 0.69). For both conditions, the values increased continuously over the sampling period (except for an outlier with small clams stranded in August 2008), from about 13 to 20 mm. Higher fluctuations in size were observed toward the end of the period, when densities were among the lowest.

Relationship to waves
The density of dying clams showed negative relationships to wave height and power for all lag periods, from 0 to 4 days, and no apparent relationship to wave period for time lag. Indeed, based on the model selection process, the significant factor was the wave height on day lag -3, i.e., the wave heights three days before the samplings. Applying the upper boundary selection resulted in a decreasing exponential relationship between the number of dying individuals and increasing wave height (Fig 3). The selected Generalized Additive Model (GAM; F = 12.07, effective d.f. = 3.73, n = 19, p<0.001) explained 80.2% of the data deviance. The size of stranded individuals showed no relationship to any physical parameter (global F 6,162 = 1.535, p = 0.170).

Relationship to the living population structure
A sharp decrease of about 15 times in the density of the T. mactroides population was recorded by Turra et al. [32] between summer and winter 2007 (Fig 2C). The shell length increased continuously from late summer, until reaching its highest mean values in winter 2008, i.e., an inverse pattern to the density (Fig 2D). In the month just prior to the period when the stranded individuals were sampled (August 2007), the buried living individuals had a mean shell length of 13.0 mm; in August 2008, their mean length was 20.0 mm.
The size of stranded individuals showed a strong relationship to the size of the living population (F 1,5 = 43.47, p<0.01, r 2 = 0.90), and the difference in the mean lengths was not significant (paired t test; t 6 = -1.74, p = 0.11). The density of stranded individuals, in turn, showed no relationship to the population density (F 1,5 = 0.01, p = 0.92).

Discussion
The fine-scale temporal monitoring employed in this study over 15 successive months, and the distinction regarding the condition of the stranded clams (dying or dead), allowed an unprecedented understanding of the frequency, magnitude and factors affecting the mortality events of beach clams, using T. mactroides in Caraguatatuba Bay as a study case.
The overall characterization of the dying and dead stranded individuals strongly supported the interpretation that stranded individuals derive from the intertidal/shallow-subtidal population stratum. This evidence includes the simultaneous presence of dying and dead stranded individuals in the upper intertidal zone, and that most of the dead individuals had their valves closed. In addition, the stranded individuals belonged to the same cohort as the buried intertidal individuals, given the close relationship of their size structure, since individuals from subtidal areas are on average much smaller than those from the intertidal zone [33]. This result highlights the ecological and behavioral components of this process, such as intertidal thermal stress, the swash dynamic, and the implications of passive displacement and relocation in this environment.
Distinguishing between dying and dead individuals showed that there is a significant relationship between the two categories. The density of dying individuals corresponded better to environmental factors, probably because they die over a short time period, which doubtless presents pseudoreplication by the sampling method used. However, separation by size structure must be done carefully. The mean dimensions of shells of dying clams were consistently larger than those of dead clams, reflecting the greater ability of larger individuals to tolerate exposure to air, since they are more resistant to desiccation [40]. Therefore, working only with dying clams may introduce biases into the size-structure information. The temporal patterns indicated that stranding of T. mactroides in Caraguatatuba Bay occurs frequently but irregularly. Even though the intensity of these events decreased considerably near the end of the study, stranding occurred in nearly the entire monitoring period, with large numbers of stranded individuals appearing sporadically. Mean daily numbers of stranded clams during October and November 2007 reached an impressive~150 ind.m -2 . The stranding events were widely distributed along this 4-km stretch of the bay, as indicated by the low standard errors of most density values.
Although quantitative data for mass-mortality events of sandy-beach species are practically nonexistent, the mean density of stranded individuals of T. mactroides, 27.09 ind.m -2 , seems high. This is particularly true considering the frequency and spatial extent of these events in Caraguatatuba Bay during the study period, but also comparing these values to estimates for the T. mactroides population in this same bay during a previous period (2003)(2004), when the population density was about 14 ind.m -2 [33]. The overall mean density of stranded individuals of T. mactroides was higher than densities recorded for Tivela isabelleana (4-15 ind.m -2 ) during mortality events, but lower than recorded for the yellow clam Mesodesma mactroides (over 60 ind.m -2 ), both on Argentinean beaches [3,13]. The broad across-and along-shore distribution of T. mactroides in the extensive Caraguatatuba Bay (~16 km long) indicates that a large amount of carrion derived from the stranded individuals becomes available to the beach-surf/ swash system and the adjacent terrestrial environment, since these clams are eaten by several marine species and the stranded individuals are also likely to attract a wide range of scavengers to the beach [20,32,50].
Burrowing clams need a wet substrate to allow the foot to penetrate and anchor before pulling the shell down, and also for the shell to remain in the same location until burrowing is completed, avoiding rolling by waves (sensu [39]). This situation is more critical for T. mactroides, which in contrast to some bivalve species that burrow quickly (e.g. Donax) [40], needs considerable time to bury and escape wave action. When burrowing in the upper swash zone of Caraguatatuba Bay, T. mactroides constricts the valves, squirts water out every 4 seconds, and buries completely only after about 20 to 30 squirts (A. Turra, pers. observation). Considering the difficulty of burying as well as the effect of low-energy waves in adding sediment to the beach [40], presumably, individuals digging out during low-energy events are more likely to be carried higher up on the beach. The surviving clams that reach the upper beach zone are then unable to return to the intertidal wet zone, because of the low wave energy that reaches that portion of the beach (weak backwash). In addition, the dry compact sand in the upper intertidal zone precludes burrowing and survival of these stranded individuals. In short, clams that are unable to bury deep enough in time would be washed up to the upper intertidal zone, by the action of an accretive wave regime, where they strand and desiccate.
This process accounts well for the negative relationship between mortality and wave height recorded here: smaller waves favor the deposition of individuals higher on the beach. Meanwhile, the higher the wave the farther from the beach it breaks, and the stronger the backwash [51], which may combine to prevent these burying intertidal individuals from being stranded. The negative effect from the harsher swash climate toward reflective conditions (lower waves), as postulated by the Habitat Harshness Hypothesis [52,53], has been widely documented for intertidal macrofauna, as also observed for this T. mactroides population, which is less dense in reflective than in dissipative conditions [33]. The negative relationship between mortality and wave height identified here is notable because stranding of mollusks, mainly bivalves, is commonly attributed to storms [3,6,18,40]. In fact, extreme conditions are prone to result in high mortality of a variety of organisms, but the present observations show that as a rule, higher mortality events were actually associated with previous low-wave-energy events. This highlights the importance of studies on finer time scales to understand certain important ecological relationships, especially in dynamic environments, and the need for these assessments in order to improve understanding and management of coastal resources.
There was no significant relationship between mortality intensity and population density, but the density of the living population will, almost as a matter of course, influence the density of stranded individuals. Data on the living population were not available with the same periodicity as the stranding data, and this lack may have resulted in the failure to identify this effect. Still there was some evidence of this relationship, most clearly in the sharp decrease in the density of stranded individuals toward the end of the study period, which was associated with a decrease in density of the living population.
For the same species and area, Turra et al. [32] found a similar natural mortality index (Z) for periods of very different densities, and showed that density-dependent processes were important in governing population dynamics. However, with respect to stranding density, the largest stranding events coincided with high population densities. The analysis of the size structure supported this suggestion, showing that the size of stranded individuals primarily reflected the size structure of the intertidal population stratum. However, whatever the population density, other factors will determine whether and which fraction of the population will die. For example, in the present study, the three cases of highest mortality (Fig 3) occurred in October and November 2007, and so may have been related to a period of high population densities. Nevertheless, the lowest mortality events were related to greater wave heights, even when the population density was quite high. Thus, the influence of wave height was very consistent, and it is reasonable to hypothesize that varying population densities would result, in the quantile regression analysis (or "envelope"), in different values for the intercept and slope of the relationship between wave height and stranding density.
Many other factors may be related to mass mortality in mollusks. A fairly well-recognized factor is the effect of parasitism, which may be of high or low magnitude, causing the death of a small portion or the collapse of a population [10,13,54]. Fiori et al. [13] suggested that parasitism by a protozoan may have caused the mass-mortality events of the yellow clam Mesodesma mactroides along its entire distribution range. This hypothesis was based on a mortality event in an Argentinean yellow-clam population infested by a coccid protozoan, which caused necrosis in the gills and stomach of some infested individuals (21.4%). Denadai et al. [55] observed parasitism in gonads of T. mactroides in Caraguatatuba Bay in the 2003-2004 period. Trematodes were found in about 10% of the individuals, and caused total infertility in some of them. In contrast to the episodic mortality events that dramatically decrease populations of M. mactroides [13], the mortality events of T. mactroides were more continuous in Caraguatatuba Bay. Thus, given the hydrodynamic and population-density aspects observed here, parasitism must not have been a main factor. However, a more-specific evaluation of parasites and pathological bacteria in this population in future mortality events is recommended, since it is a recognizable factor affecting mortality. Also, over a wider spatial scale there might be other environmental factors, such as beach slope and sediment texture, that may influence the stranding and mass mortality of T. mactroides, and could be important to understand and manage these events.
It is important to understand these high-mortality events, since large numbers of stranded clams may cause socioeconomic and environmental problems for the local residents. Harvesting of this species decreased by one-third during high-mortality periods [34]. The strong unpleasant odor of the dead clams accumulated on the sand affects tourism and raises publichealth concerns, because the clams are decomposing along the most intensively used part of the beach. Strategies to utilize this very high biomass as a solution for this environmental and socioeconomic problem should be discussed and implemented. The use of meat and shells in animal feeds for aquaculture, biodiesel production from animal fat, and use of calcium carbonate for agriculture or as an aquarium substrate are some alternatives. Nevertheless, it is necessary to undertake long-term monitoring of some of these populations to better comprehend, and predict, the dynamics of density variation.
In summary, stranding and mass-mortality events of T. mactroides were frequent and widespread along the extensive beaches of Caraguatatuba Bay in 2007-2008. The size structure of the stranded clams supports the supposition that the death of stranded individuals occurs in the intertidal and shallow subtidal zone, highlighting the importance of ecological and behavioral components to mortality processes. The wave height three days before the samplings explained much of the variation in mortality events of T. mactroides; the higher the waves, the smaller the number of stranded individuals. This is likely due to the accretive nature of lowenergy waves, so that individuals in the process of digging out would be easily carried upward, while higher waves break seaward and have a stronger backwash, preventing the clams from stranding on the uppermost part of the beach. This is an important finding for the ecology of sandy-beach bivalves, and future, temporally refined studies should examine this phenomenon further. Similarly, long-term studies on the mortality and living population of clams, concurrently with environmental features, are important for appropriate management of this resource.