Comparison of electronic length frequency analysis (ELEFAN) for estimation of growth parameters for lollyfish, Holothuria (Holothuria) atra and sand sea cucumber, Holothuria (Thymiosycia) arenicola (Holothuroidea: Echinodermata) in the north Arabian Sea, Pakistan

To estimate the growth parameters, the non-seasonal von Bertalanffy and Hoenig seasonal von Bertalanffy models were fitted to the length frequency data of Holothuria atra and H. arenicola obtained from coastal areas of Karachi in the northern Arabian Sea, Pakistan between January and December 2014. The Hoenig seasonal von Bertalanffy growth parameters were estimated as L∞=36.1 cm total length (TL), K=0.75 year-1 for H. atra and as L∞=34.9 cm TL, K=0.70 year-1 for H. arenicola. H. atra individuals reached 73.9% of their maximum total length at the one year old. For H. arenicola it was calculated as 72.7%. Monthly mean growth rate of H. atra at the same ages was calculated higher than H. arenicola. The seasonal oscillation in growth rate for H. atra (C=0.50) was larger than it was for H. arenicola (C=0.37). The time of the year when the growth is slowest corresponded to the middle of April (WP=0.30) in H. atra and at the beginning of August (WP=0.60) in H. arenicola may be the result of the extended both reproduction and poor nutrition periods due to monsoonal effects on the marine environment. The relatively high growth rates (K>0.7 year-1) of these two species may have important implications for high survival rate, particularly in environmental condition where cause biological stress and marine confusion but may also increase their potential as a candidate species for aquaculture. Please cite this paper as follows: Ahmed, Q., Bilgin, S., Ali, Q. M. (2020). Comparison of electronic length frequency analysis (ELEFAN) for estimation of growth parameters for lollyfish, Holothuria (Holothuria) atra and sand sea cucumber, Holothuria (Thymiosycia) arenicola (Holothuroidea: Echinodermata) in the north Arabian Sea, Pakistan. Marine Science and Technology Bulletin, 9(2): 145-158. * Corresponding author E-mail address: sbrbilgin@hotmail.com (S. Bilgin) Ahmed et al. (2020) Marine Science and Technology Bulletin 9(2): 145-158


Introduction
A total of 25 sea cucumber species have been recorded from the coastal waters of Pakistan, so far and seven of them are documented as commercially important species (Tahera and Tirmizi, 1995;Purcell, 2010;Purcell et al., 2012;Ahmed andAli, 2014, 2020;Ahmed et al., 2016aAhmed et al., , 2020bMoazzam and Moazzam, 2020). Currently, sea cucumber fishing occurs all over the world, but, there is no traditional sea cucumber fisheries in Pakistan (Siddique and Ayub, 2015;Moazzam and Moazzam, 2020). Two of the commercially important sea cucumber species are lollyfish, Holothuria (Holothuria) atra Jaeger, 1833 and sand sea cucumber, Holothuria (Thymiosycia) arenicola Semper, 1868 found in the coastal waters of Pakistan.
H. atra known also as the black sea cucumber is found in the tropical Indo-Pacific region, its range extends from the Red Sea and East Africa to Australia. It inhabits the inner and outer flats, back reefs, shallow lagoons, sand-mud and rubble, and sea grass beds between 0 and 20 m. In few Pacific Island nations, the body wall, intestines and/or gonads are consumed in traditional diets or in times of hardship (Purcell et al., 2012).
H. arenicola is also found from the Indo-Pacific to the tropical Western Atlantic. The species is reported from tropical Australia, Mozambique, the Red Sea, Ascension Island, Virgin Islands, Antigua, Barbados, Tobago, Belize, Bermuda, Tortugas, Jamaica, Puerto Rico, U.S. (Mosher, 1980;WoRMS, 2020). These two holothurians are found on the seabed, in shallow waters on reefs and under intertidal rocks, sand flats and in seagrass meadows at depths of up to 20-30 m (Mosher, 1980;Pourvali et al., 2014;WoRMS, 2020). H. atra stocks in many regions are now being exploited and are often marketed under the name bêche-de-mer or trepang which means "cooked and dried sea cucumbers" (Purcell et al., 2012).
Relatively little information is available on the biology and ecology of these sea cucumbers belong to the class Holothuroidea in the northern Arabian Sea coasts of Pakistan. There are some works about the growth of sea cucumbers; including the weight-length relationships (WLRs) and condition factor (CF) based growth features in holothurian species such as Ohshimella ehrenbergii Selenka, 1867, H. arenicola, H. atra, H. pardalis Selenka, 1867and H. verrucosa Selenka, 1867 from this region (Siddique et al., 2015;Ahmed et al., 2018aAhmed et al., , 2018b. In addition to these, a detailed study on the population dynamics of sea cucumbers has been carried out on the H. arenicola stocks in Manora and Buleji rocky shores in the northern Arabian Sea, Pakistan so far (Siddique and Ayub, 2015). It was reported that an increase of the gonad index (GI) of H. arenicola was observed during spring and early summer, followed by a decrease in autumn and winter, which showed the spawning followed by resting phase and also the GI values were reported as a significant negative correlation with salinity and non-significant correlation with temperature (Siddique and Ayub, 2015). Moreover, unimodal distribution, bimodal distribution and three modes were determined for H. arenicola by seasons and areas (Siddique and Ayub, 2015). In other geographical regions, previous studies have been conducted for H. atra to determine reproduction biology features including sexual and asexual reproduction from Heron Reef (Harriott, 1982(Harriott, , 1985, spawning season, size at sexual maturity and fecundity from Red Sea (Chao et al., 1994;Abdel-Razek, 2005), ecology including distribution, habitat preference and feeding activity from Rongelap Atoll (Bonham and Held, 1963), Bahamas (Mosher, 1980), and Gulf of Mannar (Asha et al., 2015) and also environmental variability, population biology and fishing pressure from South Tarawa Lagoon, Republic of Kiribati (Tamaroa, 2010).
Age determination of holothuroid echinoderms is difficult because of their soft-bodied, polymorphic shape, and the small size of hard structures that could exhibit growth rings (Watanabe et al., 2014;Sun et al., 2019). Analyses of length frequency data (LFDA) can be used for estimating age classes and growth parameters (Pauly and David, 1981), and can also be applied to calculate growth of holothuroid echinoderms (Olaya-Restrepo et al., 2018). An additional challenge for estimating growth of holothuroid echinoderms arises due to the seasonal growth pattern, and a modified von Bertalanffy growth model has been developed to incorporate this (Hoenig and Hanumara, 1982). Seasonal growth patterns have been reported for different Holothuria species, such as H. arguinensis Koehler andVaney, 1906 (Olaya-Restrepo et al., 2018), Isostichopus badionotus (Selenka, 1867) (Poot-Salazar et al., 2015 and Australostichopus mollis (Hutton, 1872) (Morgan, 2012). However, no seasonal growth information exists for the mentioned unexploited Holothuria species inhabiting the Arabian Sea. The aim of this study was to investigate seasonal growth rate of H. atra and H. arenicola in the north Arabian Sea.

Study Area and Sampling
Holothuria atra (n=221) and Holothuria arenicola (n=258) specimens were collected on monthly basis from Buleji (24°50'20.41" N, 66°49'24.15" E) and Sunehri (24°52'33.49" N, 66°40'40.20" E) ( Figure 1) coasts from January to December 2014 except June and July due to rough sea in monsoon season. The sample was collected from intertidal zone by hand-picking through forceps at low tide. Collected specimens were kept alive in water filled containers and then were transported to the laboratory and shifted in well aerated aquaria. For taxonomic studies and identification, morphological features were examined and microscopic studies were conducted. Identification of species of sea cucumbers was based on characteristics body morphology of the sea cucumbers, using the guidebooks from Selenka (1867) and Samyn et al. (2006).
Ossicles were taken from three positions (dorsal and ventral body walls, and tentacles); wet mounts were prepared by placing a small piece of skin tissue on slide and adding few drops of 3.5% bleach, the slides were then rinsed with drops of distilled water. The slides were examined under the Nikon LABOPHOT-2 microscope at 10×10 magnifications. Microphotography was also performed through Fujifilm 16 MP digital camera (see Ahmed et al., 2018aAhmed et al., , 2018b for more details). Length (cm) data were collected for each sea cucumber after allowing the sea cucumber to relax in water for 5 min. Total length from mouth to anus was measured with the flexible ruler.

Von Bertalanffy Growth Function Parameter Estimation
Growth in length has been described using the von Bertalanffy (1938) growth function, based on either observed or back calculated length at ages. The length frequency distribution analysis (LFDA) package is also a PC based computer package for estimating growth parameters from length frequency distributions. Version 5.0 of LFDA includes methods for estimating the parameters of both non-seasonal and Hoenig seasonal versions of the von Bertalanffy growth curve (Kirkwood et al., 2003).
The standard or non-seasonal von Bertalanffy (1938) growth function (VBGF) is: Seasonal growth or five parameters von Bertalanffy growth model (5 Parameters VBGF) was described using the Somers's (1988) version of the VBG equation: where, Lt is length at age t, L∞ is the asymptotic length to which the sea cucumber growth, K is the growth-rate parameter, t0 is the nominal age at which the length is zero, C is the relative amplitude (0 < C < 1) of the seasonal oscillations, tS is the phase of the seasonal oscillations (-0.5 < tS < 0.5) denoting the time of year corresponding to the start of the convex segment of sinusoidal oscillation.
The time of the year when the growth rate is slowest, known as the winter point (WP), was calculated as: Seasonal VBG curves were fitted to the length distributions after first indicating a range of values of L∞ and K and reducing iteratively the range to maximize the goodness of fit (Rn) of the curves to the data. Rn was calculated as: where ASP is the available sum of peaks, computed by adding the best values of the available peaks, and ESP is the explained sum of peaks, computed by summing all the peaks and troughs hit by the VBGF curve. Analysis of the length data were fitted to length frequency distributions grouped in 2 cm total length size classes using the Electronic LEngth Frequency ANalysis (ELEFAN) procedure in the PC-based computer package Version 5.0 of Length-Frequency Distribution Analysis (Kirkwood et al., 2003). The ELEFAN procedure first restructures length frequencies and then fits a VBGF curve to the restructured data. Both seasonal and non-seasonal VBGF curves were fitted to the seasonal length distribution after providing a range of values for the parameters to be estimated and then iteratively reducing the range until the goodness of fit of the curve to the data is maximized.

Reliability of Growth Parameter Estimates
Having estimated a set of growth parameters, one would like to evaluate their reliability. A possible test is the so-called phiprime test (Φ') known as growth performance index. This test is based on the discovery by Pauly and Munro (1984) that Φ' values are very similar within related taxa. So, the growth performance comparisons were made using the growth performance index (Φ') which is preferred rather than using L∞ and K individually (Pauly and Munro, 1984) and is computed as:

Seasonal and Non-Seasonal Von Bertalanffy Growth Parameters
The Hoenig seasonal and non-seasonal VBGF curves parameters obtained from the LFDA for each species are summarized in Table 1. The Rn value of the non-seasonal growth curve for both H. atra and H. arenicola improved when the seasonal growth curves were fitted ( Table 1), suggesting that, at least for our data, these two sea cucumbers exhibit seasonal growth patterns. This seasonality in growth was also apparent in the results of the relative amplitude values (C=0.50 for H. atra and C=0.37 for H. arenicola) of the seasonal oscillations and in Figure 4B and Figure 5B where sinusoidal pattern could be observed in the Hoenig seasonal von Bertalanffy growth curve.
The slowest growth period started at the middle of April for H. atra (WP=0.30; Figure 6). For H. arenicola, however, the start of slow growth period was at the beginning of August (WP=0.60; Figure 6).
The calculated growth performance indices (Φ'; This fast growth characteristic of small individuals was also apparent in the growth curves in Figure 4 and Figure 5, where it could be observed that the slight slope in the larger individuals compared to smaller sea cucumber leading to small individuals grew faster than large ones. Monthly mean growth rate of H. atra was calculated as 1.18 cm/month for 0-365 days (0 age class), 0.97 cm/month 365-730 days (1 age class), 0.45 cm/month for 730-1095 days (2 age class) by using estimated average total length values in the ages.

Discussion
To the best of our knowledge this is the first study to calculate the Hoenig seasonal growth curve parameters by using ELEFAN, to length-frequency data and age-length key of two holothurians, H. arenicola and H. atra. Our LFDA results showed that both H. arenicola and H. atra have similar lengthfrequency distribution and modes. Different spatial and temporal variation of LFDA and/or modes was reported for H. atra from Manora and Buleji in the previous study (Siddique and Ayub, 2005). According to Siddique and Ayub (2005), the small-sized individuals (<12 cm) were equally abundant at both sites, the medium-sized individuals (12-20 cm) were more (37.7% of the total population) at Buleji than at Manora (20.6%) and the larger-sized individuals (>22 cm) were higher in numbers at Manora (21.2%) than at Buleji (2.4%). Moreover, bimodal distribution (consisting of small and medium-sized individuals in summer and autumn), three modes (winter and spring) from Manora and unimodal distribution from Buleji were reported. The mean total length of examined individuals in Buleji (mean length: 11.6±3.5 cm, between 3.8 and 37.6 cm) was also reported as higher than Monara (mean length: 14.1±7.6 cm, between 2.7 and 41.6 cm) during the period from April 2011 to November 2012. In our study, the mean total length (mean: 23.3±0.30, between 13.0 and 33.0 cm) of examined individuals were about 10-11 cm larger than the study of Siddique and Ayub (2005) and also, we have not observed small-sized individuals (<12 cm). In the present study, the sea cucumbers smaller than 10-12 mm could not be obtained in the sampling area both two species, most probably due to the differences in behavioral characteristics (e.g. burrowing cycle) by size. Mercier et al. (2000) reported that juveniles H. scabra ≤10 cm burrowed at sunrise and surfaced at sunset, whereas individuals >10 cm burrowed and surfaced a few hours earlier (in the middle of the afternoon). According to these results described above, it can be concluded that the growth and length composition for examined species in the research region differ significantly between years and between regions due to the marine environmental variability such as site health, e.g., water quality (dissolved oxygen and nutrients), hydrodynamics (monsoonal effect, current, water movement). Secondly, because of the different natural mortality and migration or recruitment levels. Migration of H. atra from one area to another was reported by different researchers due to rolling with the water movement (Tamaroa, 2010) and also due to lack of dissolved oxygen and nutrients in any environment (Uthicke, 2004).  The observed growth rates (K>0.70) indicate that examined two species achieve asymptotic size quickly, even faster some of the other holothurians, such as I. badionotus (K=0.2), I. fuscus (K=0.18), H. pulla (K=0.24). On the other hand, similar growth parameters were also reported for many other holothurian species such as H. arguinensis (K=0.88), I. badionotus (K=0.7) ( Table 2). The general assumptions for different length-based growth models for the growth determination of holothurians reported as (Olaya-Restrepo et al., 2018): i) species reproduce and recruit once per year; ii) species reach a specific size; and iii) species grow isometrically (b=3). Most commercial sea cucumber species fulfill the first two assumptions, but usually they do not exhibit isometric growth (Olaya-Restrepo et al., 2018). For H. atra and H. arenicola, negative allometric growth has been reported (Siddique et al., 2015;Ahmed et al., 2018b).
However, much larger species such as I. badionotus are long lived and usually have relatively low values for K. In such cases it could be expected that life span is extended considerably. For I. badionotus in Sisal, at asymptotic population size was reported ten years while K was as low as 0.20 year -1 (Poot-Salazar et al., 2015). Faster growth rates may have important and positive effects on the mortality of holothurians such as Australostichopus mollis (Hutton, 1872), primarily under high stress conditions and in disturbed environments (Morgan, 2012) and so, the fast growth could be considered as an adaptive advantage. Namely, larger sea cucumbers would be capable of moving and stabilizing themselves in more mobile sediments, such as sandy substrates, and they would also be more likely to escape predation (Mercier et al., 2000;Olaya-Restrepo et al., 2018). In our study, calculated K (0.75 year -1 ) and L∞ (36.1 cm) values for H. atra were found to be larger than the K (0.11 year -1 ) and L∞ (32.4 cm) values reported in the previous study, conducted between 1975 and 1976 on the seaward reef bench of Ananij Island, Enewetak Atoll, Marshall Islands (Elbert, 1978) (see: Table 2). This difference between our study and Elbert's study, which was carried out about 42 years ago, may have resulted from the difference between biotic and abiotic factors such as prey availability, predators, genetic variation, salinity, habitat structure, etc. and also may be due to the difference of the method used to calculate of the growth. When there is a seasonal growth pattern for holothurian species belonging to the same family in a geographical region, the estimations of L∞ and K may differ significantly between the seasonal and nonseasonal models. In the present study, the non-seasonal VBG model provided realistic results. However, when seasonality was included (with the Hoenig model), more reliable values were obtained, which confirmed the seasonality in the growth of H. arenicola and H. atra.  The Rn and C values with visual growth curves both H. arenicola ( Figure 4B) and H. atra ( Figure 5B) evidenced that these two species exhibited marked seasonality in growth.
Seasonal growth pattern was also reported from different  (Chao et al., 1994), New Caledonian lagoon (Conand, 1993), and at Heron Reff, Greak Barrier Reff (Harriott, 1985). Temperature, salinity, food abundance and photoperiod have all been known as factors regulating holothurian reproductive cycle. The coast of Pakistan is for most of the year influenced by high-salinity surface water (36-38‰) and the sea surface temperature during summer is 28-30°C while during winter, it is 21-24°C (Siddique and Ayub, 2015). Furthermore, there is monsoonal variations in nutrient concentration due to monsoon winds and fluctuations of rainfall (Hussaina et al., 2010). Such variations may also be related to productivity and avai1ability of food and to the reproductive cycle and growth of H. arenicola and H. atra and also other sea cucumber species in the north Arabian Sea, Pakistan (Ahmed et al., 2020a). All these could be considered as the factors forcing seasonality in growth both H. arenicola ( Figure 4B) and H. atra ( Figure 5B) in the study area. The variations in the period of the slowest growth time (WP) generally related to environmental condition which cause biological stress, physiological conditions of the marine animal, fullness of stomach and gonads stages (Bilgin et al., 2009a(Bilgin et al., , 2009bAhmed et. al., 2016bAhmed et. al., , 2020a. Temperature also plays an important role in the reproductive events and the abundance of food, and therefore especially the water temperature indirectly effective on the WP of the marine organisms. The seasonal and/or monsoonal variations in the mean condition factors of H. arenicola, H. atra, H. pardalis and H. verrucosa was reported from the coasts of Karachi, Pakistan (Siddique et al., 2014;Ahmed et al., 2018aAhmed et al., , 2018b and also the gonad index (GI) based reproduction time for H. arenicola was reported during spring and early summer from these regions (Siddique and Ayub, 2015). The period of slowest growth for H. arenicola corresponded to June when the highest GI value occurs (see Siddique and Ayub, 2015). The slowest growth period for H. atra, however, was estimated to be April, which is the period with a relatively low water temperature and with the most probably when reproductive cells (gonad development) begin to mature or develop (Harriott, 1985;Conand, 1993;Chao et al., 1994;Abdel-Razek et al., 2005). Since the growth rate of holotorians depends on the food availability and quality of the habitat and/or organic matter content of sediment, dissolved oxygen concentration in the water column and in sediment (Tamaroa, 2010) and also reproduction events, the slow growth of H. atra in April may be the result of the extended both poor nutrition and reproduction periods growth but this was not tested in this study. Therefore, the effects of these factors on growth should be examined in the future studies.

Conclusion
Together with this study seasonal growth was first pronounced for both H. arenicola and H. atra. Although the causes of this seasonality are not tested in this study, the effect of fluctuations in the amount of nutrient and maturity levels are probably important for this situation. The effects of the climatic events as monsoon rain and reproduction activities on growth should be studied in detail in order to provide data for holotorian fisheries management in the north Arabian Sea, Pakistan.