Identification of fish families and species from the western Arabian Gulf by otolith shape analysis and factors affecting the identification process

Additional keywords: linear discriminant analysis, sampling procedure, species identification, wavelet transformation.

Otolith shape analysis is a well-established method of delineating fish stocks, populations, species and higher taxonomic groups. It is fairly cost-effective and requires only otolith images containing information about shape, outline and landmarks. The analysis programs are usually built in a free software environment with introductory protocols (e.g. the shapeR package, see https://cran.r-project.org/package=shapeR, in R, ver. 3.4.2, R Foundation for Statistical Computing, Vienna, Austria, see https://www.R-project.org; Libungan and Pálsson 2015), which facilitates wide application. The requirement for digital images only further enables large electronic catalogues to be established and automated taxon identification systems to be developed that can be accessed by the general public (e.g. the Anàlisi de Formes d’Otòlits (AFORO) database, http://isis.cmima.csic.es/aforo/, accessed 31 July 2018; Lombarte et al. 2006).

Although otolith shape analysis is a promising tool for the identification of fish groups, several factors may hinder its successful use, the first being the selection of available analytical methods. Differences in otolith shape have been used to identify teleost species for over 130 years (Nolf 2013; Lombarte et al. 2018) and numerous analytical approaches have been developed (e.g. Tuset et al. 2003, 2016; Ponton 2006; Wakefield et al. 2014; Libungan and Pálsson 2015). The performance of these methods in identifying fish groups is often case specific (Ponton 2006), and it remains unclear whether incorporating different sources of information would enhance the identification performance of these methods. Second, identification performance may be affected by the number of candidate species. Lombarte et al. (2018) suggested an inverse relationship between the correct identification rate and the number of candidate species based on a summary of six studies. However, these studies differed in methodology and collection sites, and the relationship has not been tested in a systematic manner using the same methodology in a single study. Third, the sample sizes used for model construction vary considerably among studies, from 10 otoliths for each species in the study of Salimi et al. (2016) up to 136 otoliths for the goby Sicyopterus lagocephalus in the study of Lord et al. (2012). Fourth, otolith shape may change allometrically with somatic growth (Monteiro et al. 2005; Vignon 2012). Therefore, sufficient sample sizes with proper sampling procedures, such as size-stratified random sampling, could ensure that variations in otolith shape due to ontogenetic changes are included in the samples. It is of practical importance to know how different sample sizes and sampling procedures may affect identification performance and to design a best practice that balances the trade-off between sampling costs and identification precision.

The primary objective of this study was to identify fish families and species from the western Arabian Gulf using otolith shape analysis. Second, we tested the effects of several factors on identification performance, including the use of different data sources, sex as a possible confounding factor (Cardinale et al. 2004), the number of candidate species, different sample sizes and different sampling procedures with regard to somatic sizes. Four data sources were extracted from otolith shape: otolith measurement, shape indices (Tuset et al. 2003), coefficients from wavelet transformation (Libungan and Pálsson 2015) and relative otolith length. Fisher’s linear discriminant analysis was used to analyse the multivariate data and the classification of fish families and species.

Fish specimens were collected from bottom trawl surveys (Rabaoui et al. 2015, 2017) and visits to major landing sites under a landing site monitoring program (Rabaoui et al. 2017). The trawl surveys involved a commercial outrigger trawler in the western Arabian Gulf in 2013 and 2016. The landing site monitoring program examined 12 fish and invertebrate stocks of major commercial importance in 2013 and 2014.

In this study, we selected 16 species from 5 families of commercial interest: emperor fish Lethrinus nebulosus, Lethrinus lentjan, Lethrinus microdon and Lethrinus borbonicus from Lethrinidae, threadfin bream Nemipterus bipunctatus, Nemipterus japonicus, Nemipterus peronii and Nemipterus randalli from Nemipteridae, groupers Epinephelus coioides, Epinephelus areolatus and Cephalopholis hemistiktos from Serranidae, seabream Argyrops spinifer, Rhabdosargus haffara and Sparidentex hasta from Sparidae and lizardfish Saurida macrolepis and Saurida tumbil from Synodontidae. After thawing, the total length and weight were measured to the nearest 1 mm and 0.1 g. The largest pair of otoliths, the sagittae, was removed from the fish, cleaned in fresh water, dried in air and stored in sealed plastic vials.

A stratified random sampling scheme was used to fully represent the size distribution of the population. The fish were ranked by somatic size and classified into three size strata, namely small, medium and large, covering the lower, middle and upper thirds of the somatic length rank respectively. Between 17 and 60 specimens were randomly drawn from each stratum. In all, 1953 fish specimens were collected, with sufficient coverage in length that fish were sampled from juveniles to large adults (Table 1).

Table 1.  Summary table for the fish collected from the western Arabian Gulf
L_T, total length

The otoliths of some fish species become greatly curved towards the external face and their three-dimensional (3D) structures become more apparent as the fish grows. Thus, it was not feasible to obtain otolith photographs in which all the otoliths were positioned in a perfect horizontal plane with the sulcus acusticus facing the observer, as done by Lombarte et al. (2018). Alternatively, the otoliths were placed horizontally with the external plane of the otolith facing towards the observer (Fig. 1), as described by Libungan and Pálsson (2015). The otoliths were observed against a black background with a reflected light source under a dissection microscope attached to a digital camera (DP72; Olympus, Center Valley, PA, USA). The images were captured and digitised using CellSens Standard software (ver. 1.5, Olympus).

Fig. 1.  Otolith photographs for 16 fish species from 5 families: groupers Epinephelus coioides, E. areolatus and Cephalopholis hemistiktos (otolith lengths 15.4, 8.8 and 8.4 mm respectively) from Serranidae; seabream Argyrops spinifer, Rhabdosargus haffara and Sparidentex hasta (otolith lengths 17.5, 7.9 and 11.4 mm respectively) from Sparidae; lizardfish Saurida macrolepis and S. tumbil (otolith lengths 5.9 and 12.8 mm respectively) from Synodontidae; emperor fish Lethrinus nebulosus, L. lentjan, L. microdon and L. borbonicus (otolith lengths 12.9, 13.2, 7.3 and 5.6 mm respectively) from Lethrinidae; and threadfin breams Nemipterus bipunctatus, N. japonicus, N. peronii and N. randalli (otolith lengths 7.0, 10.2, 6.9 and 6.0 mm respectively) from Nemipteridae.

Extraction of the otolith outline, contour smoothing, measurement, generation of shape coefficients (wavelet transformation and normalised elliptic Fourier transformation) and standardisation of the shape coefficients for fish length were completed following the protocol of Libungan and Pálsson (2015) using the shapeR package (Libungan and Pálsson 2015) in R (ver. 3.4.2, R Foundation for Statistical Computing). Wavelet transformation resulted in a better approximation of the otolith outline (0.4% deviation from observed outline at Level 5) than did normalised elliptic Fourier transformation (1.8% at the 25th harmonic), and only the standardised shape coefficients from wavelet transformation were used in the classification of fish families and species.

Four measurements were available after running shapeR, namely the otolith length (L_O; mm), otolith width (W_O; mm), otolith perimeter (P_O; mm) and otolith area (A_O; mm²). Then, five dimensionless shape indices were calculated according to Tuset et al. (2003) to provide an alternative source for otolith shape, namely form factor (4πA_O ÷ P_O²), roundness (4A_O ÷ (πL_O)²), circularity (P_O² ÷ A_O), rectangularity (A_O ÷ (L_O × W_O) and ellipticity ((L_O – W_O) ÷ (L_O + W_O)). Further, relative otolith length (L_O,Rel) was calculated as L_O ÷ L_T, where L_T is total somatic length (mm).

Fish families and species were classified using Fisher’s linear discriminant analysis (LDA) applied to different datasets: 4 otolith measurements (otolith length, otolith width, otolith perimeter and otolith area), 5 shape indices (form factor, roundness, circularity, rectangularity and ellipticity), 26 wavelet coefficients and relative otolith length. Different LDA models were constructed using five combinations of available datasets: (1) otolith measurement; (2) shape indices; (3) wavelet coefficients; (4) shape indices plus wavelet coefficients; and (5) shape indices plus wavelet coefficients plus relative otolith length. Shape indices were calculated from otolith measurements; thus, these two datasets were of the same source and were not combined. Two situations were examined: (1) the families and species are unknown and the objective is to identify the families or species; (2) the family is known a priori and the objective is to identify the species given the family. The performance of LDA in identifying families and species was evaluated by the correct prediction rate estimated by jack-knife (leave-one-out) method. The LDA was run using the R package MASS (see https://cran.r-project.org/package=MASS; Venables and Ripley 2002).

The effects of somatic length and sex on otolith shape were evaluated using the permutation test for canonical analysis of principal coordinates (Anderson and Willis 2003) with 9999 permutations using the function anova.cca in the R package vegan (ver. 2.5–2, J. Oksanen, F. G. Blanchet, M. Friendly, R. Kindt, P. Legendre, D. McGlinn, P. R. Minchin, R. B. O’Hara, G. L. Simpson, P. Solymos, M. Henry, H. Stevens, E. Szoecs and H. Wagner, see https://CRAN.R-project.org/package=vegan, accessed 31 July 2018). Shape indices, wavelet coefficients and relative otolith length were the response variables, and somatic length and sex were covariates. Because E. coioides, E. areolatus, C. hemistiktos and S. hasta are hermaphrodites (the former three being protogynous and the last being protandrous; Frisch 2004) and sex identification for N. bipunctatus, N. japonicus, N. peronii and N. randalli was not feasible because of the lack of well-developed gonads, the permutation tests were applied to L. nebulosus, L. lentjan, L. microdon, L. borbonicus, A. spinifer, R. haffara, Saurida macrolepis and S. tumbil. Sex was categorised as female or male, and sex-undifferentiated individuals of these species were excluded from the analysis.

Simulation was used to evaluate the effects of different numbers of candidate species, sample size and sampling procedures on the identification of fish families and species using otolith shape analysis.

To evaluate how the number of candidate species affected the performance of otolith shape analysis in species identification, 2–14 species of the 16 species examined were randomly selected, yielding a total of 120–12 870 combinations. The LDA was conducted over all combinations using five combinations of available datasets. The correct prediction rate (%) in the null model, in which the species was correctly identified only by chance was calculated by dividing 100 by the number of candidate species.

To examine the effects of sample size and sampling procedure, the following scenarios were created and simulated:

1. Reference: the LDA was conducted based on the full sample

2. Simple random sampling: 6, 12, 15, 18, 24 or 30 individuals in total from each species were randomly selected without replacement from the sample for use in the LDA

3. Size-stratified random sampling: for each species, individuals were classified into three size classes (small, medium and large) based on total length, after which 2, 4, 5, 6, 8 or 10 individuals in total were randomly selected without replacement from each size category for each species, leading to total sample sizes of 6, 12, 15, 18, 24 and 30 for the LDA

4. Size-biased random sampling: three cases were contained in this scenario (small-, medium- and large-biased), representing three sources of bias with regard to somatic size; for each species, 15 individuals were randomly selected without replacement from the small, medium and large size categories for use in the LDA.

In all, 5000 simulations were run for each of the four scenarios described above. Selection without replacement by species and size category was performed using the ‘sampling’ package in R (ver. 2.8, Y. Tillé and A. Matei, see https://CRAN.R-project.org/package=sampling, accessed 31 July 2018).

Performance identifying 5 families and 16 species using different datasets, including 4 otolith measurements, 5 shape indices, 26 wavelet coefficients and the relative otolith length, is summarised in Table 2. Generally, no single dataset outperformed any other, and the incorporation of all available datasets (i.e. shape indices, S, plus wavelet coefficients, W, plus L_O,Rel; Table 2) resulted in the highest correct prediction rate for the identification of families and species. This combination was able to distinguish Synodontidae, Serranidae and Nemipteridae with a correct prediction rate of 89–100%. However, it was more difficult to identify Sparidae and Lethrinidae, with a prediction rate at best of 63–78%. For species identification, the yellowfin hind C. hemistiktos was generally identified with a good precision of 71–93%, whereas it was most difficult to identify the Sobaity seabream S. hasta (low precision of 6–62%).

Table 2.  Correct prediction rate for the identification of fish families and species using four otolith measurements (M; i.e. otolith length, width, perimeter and area), five dimensionless shape indices (S), 26 wavelet coefficients (W) and relative otolith length (L_O,Rel)

When the families were known a priori, the precision of identifying species increased substantially, and 13 species could be identified with a precision of >90% (Table 3). The incorporation of different datasets generally enhanced the identification of species, especially for L. microdon, N. japonicus, E. areolatus, S. hasta and S. tumbil.

Table 3.  Correct prediction rate for the identification fish species given the family is known a priori using four otolith measurements (M; i.e. otolith length, width, perimeter and area), five dimensionless shape indices (S), 26 wavelet coefficients (W) and relative otolith length (L_O,Rel)
Values in parentheses for Lethrinidae, Serranidae and Synodontidae show the median correct prediction rate for four, three and two randomly chosen species respectively

Somatic length significantly affected all datasets that included otolith shape (e.g. S+W+L_O,Rel; permutation test for canonical analysis of principal coordinates with 9999 permutations, all P ≤ 0.0001), whereas sex (female v. male) did not significantly affect otolith shape (all P ≥ 0.2272) for L. nebulosus, L. lentjan, L. microdon, L. borbonicus, A. spinifer, R. haffara, S. macrolepis and S. tumbil (Table 4).

Table 4.  Analysis of variance (ANOVA) table for the permutation test for canonical analysis of principal coordinates with 9999 permutations
Sex (female v. male) and somatic length were the two factors tested simultaneously, with shape indices, wavelet coefficients and relative otolith length the response variables

The performance of otolith shape analysis in identifying fish species clearly declined with an increasing number of species known (Fig. 2). LDA models from four datasets and their combinations had higher correct prediction rates than the null model. Using all sources of data led to the best performance, with the highest correct prediction rate and smallest variation (Fig. 2).

Fig. 2.  Boxplots showing correct prediction rates for fish species identification when 2–14 species are randomly chosen from a total of 16 species using (a) otolith measurements, (b) shape indices (S), (c) wavelet coefficients (W), (d) S+W and (e) S+W+relative otolith length. The grey dotted lines indicate the null mode of being correctly predicted by chance. The boxes show the interquartile range, with median values indicated by the horizontal lines; whiskers show the range. Individual symbols indicate outliers.

When identifying fish families, performance increased gradually with increasing sample size, as indicated by increasing median correct prediction rate (Fig. 3). As the sample size increased, the variation in the prediction rate declined. However, for Synodontidae, which was identified with good precision, increased samples sizes did not result in better performance (Fig. 3).

Fig. 3.  Boxplots showing correct prediction rates of identifying five families under two sampling schemes and with seven sample size levels (6, 12, 15, 18, 24, 30 and full sample size). The upper panels show results when fish were sampled by simple random sampling; the lower panels show results when fish were sampled by stratified random sampling from three size categories (small, medium and large). The boxes show the interquartile range, with median values indicated by the horizontal lines; whiskers show the range. Individual symbols indicate outliers.

When identifying fish species, increasing sample size also generally led to improved precision with less variation in the prediction rate, but this improvement was species specific (Fig. 4). For example, the correct prediction rate did not increase for S. tumbil and S. macrolepis with increasing sample size. Size-stratified random sampling led to similar or slightly better correct prediction rates, and the results were more robust to changes in sample size (Fig. 4).

Fig. 4.  Boxplots showing correct prediction rates of identifying 16 species under two sampling schemes and for seven sample size levels (6, 12, 15, 18, 24, 30 and full sample size). (a, c) Fish were sampled by simple random sampling; (b, d) fish were sampled by stratified random sampling from three size categories. The boxes show the interquartile range, with median values indicated by the horizontal lines; whiskers show the range. Individual symbols indicate outliers.

The effects of fish size-biased sampling procedures were family and species specific. When identifying families, size-biased sampling procedures produced no differences relative to random sampling and the reference case. Sparidae was the exception, where small-biased sampling had a better prediction rate than medium- and large-biased sampling procedures (Fig. 5). When identifying species, the responses to size-biased sampling differed among species (Fig. 6). The median correct prediction rate for L. borbonicus, N. peronii, C. hemistiktos and R. haffara decreased as the biases in size moved from small to large. There was a considerable decrease in the correct prediction rate for A. spinifer when the samples consisted mostly of large-sized individuals (Fig. 6).

Fig. 5.  Boxplots showing correct prediction rates of identifying five families under six sampling procedures. In the Small, Medium and Large scenarios, estimates of correct prediction rate were derived from a total of 15 samples randomly drawn from small-, medium- and large-sized individuals respectively. In the Sim.Ran and Stra.Ran scenarios, the 15 samples were selected by simple random sampling and size-stratified random sampling respectively with three size categories. In the Full scenario, 10 individuals were chosen a priori as the test sample and analyses were performed using the remaining individuals. The boxes show the interquartile range, with median values indicated by the horizontal lines; whiskers show the range. Individual symbols indicate outliers.

Fig. 6.  Boxplots showing correct prediction rates of identifying 16 species under 6 sampling procedures. In the Small, Medium and Large scenarios, estimates of correct prediction rate were derived from a total of 15 samples randomly drawn from small-, medium- and large-sized individuals respectively. In the Sim.Ran and Stra.Ran scenarios, the 15 samples were selected by simple random sampling and size-stratified random sampling respectively with three size categories. In the Full scenario, 10 individuals were chosen a priori as the test sample and analyses were performed using the remaining individuals. The boxes show the interquartile range, with median values indicated by the horizontal lines; whiskers show the range. Individual symbols indicate outliers.

By demonstrating different performances in identifying families and species among datasets, we have shown the importance of having different sources of data to fully represent the variations in otolith shape and outline. Different sources of data are not mutually exclusive, and they could be incorporated by simply adding them into the multivariate dataset for LDA. Modern approaches, such as a machine‐learning framework (Chen et al. 2018), could be promising for incorporating different datasets.

Performance identifying species declined with the number of candidate species. Therefore, incorporating information a priori to narrow down the number of candidate species could produce better correct prediction rates, because identifying species from a known family is simpler. For example, identifying species from otoliths in fossiliferous deposits, seabed sediments and archaeological remains can provide great information about historical baselines (Lin et al. 2015, 2016). In such cases, expert identification may have classified these otoliths into larger taxonomic groups, such as family or genus (Beech 2001; Girone and Nolf 2009). Otolith shape analysis can then be used to refine the taxonomic resolution.

Correct species identification relies on whether a sufficient number of species in the genus and family have been included so that the within- and between-species variations in shape are well represented in the model. Therefore, the credibility of species identification by otolith shape analysis depends strongly on the current knowledge of fish taxonomy and species delimitation. In the Arabian Gulf there are four species in the family Lethrinidae, nine species in the family Nemptereridae, which has four species in the genus Nemipterus, and four species in the family Synodontidae, which has two species in the genus Saurida (Carpenter et al. 1997). The present study has completely examined all species in the family Lethrindae and the genera Nemipterus and Saurida in the Arabian Gulf. All between-species variations in shape in these three groups have been included in our models. If these three groups are identified a priori, otolith shape analysis is a very powerful tool for identifying the species within these groups with high correct prediction rates of 83–100%.

Sex did not affect the otolith shape parameters for selected species, which is consistent with other studies (e.g. Begg et al. 2000; Cardinale et al. 2004). Strong effects from somatic length indicated that evolution in otolith shape with somatic length was not completely excluded by standardisation of otolith shape parameters for somatic size. Therefore, it is important that sampling obtain sufficiently large sample sizes covering all available size classes so as to fully represent possible variations in otolith shapes and obtain good prediction precision.

Variations in otolith shape could result from many factors, including different fish habitats (Lord et al. 2012), hearing and sound production functions (Cruz and Lombarte 2004), taxonomic and phylogenetic relationships (Lombarte and Cruz 2007; Tuset et al. 2008, p. 10; Lin and Chang 2012), and ontogenetic effects in different developmental stages (Monteiro et al. 2005; Vignon 2012). Insufficient sample sizes could lead to underestimation of the true variation in otolith shapes, and consequently could lead to reductions in correct prediction rate. Biases in the size distribution of specimens could further affect the otolith shape analysis if the species of interest exhibits strong differences in otolith shape with size, such as the king soldier bream A. spinifer. When it is difficult to obtain a large sample size, a size-stratified random sampling program could be used to spread sampling effort over a wide range of sizes to account for possible allometric changes in otolith shape.

Finally, we would like to address the issue of incomplete species lists where not all species in a given higher taxonomic level (e.g. Family) could be included in the analysis. This is especially common for recovered fossil or archaeological materials. Recent taxonomic descriptions of new species, such as the Sparidae in the Arabian Gulf (Iwatsuki 2013; Amir et al. 2014), can further increase the total number of species. Such uncertainty in the species list and the possibility of new and undescribed species may lead to biased between- and within-species variability in otolith shape, which may explain the low correct prediction rate for some sparid species, such as S. hasta. We proposed two multivariate statistical approaches for the case of incomplete species list. If the number of species included in the analysis is limited, plotting all samples over a plot after some dimension-reduction method, such as canonical analysis of principal coordinates (Libungan and Pálsson 2015) and canonical variate analysis (Lombarte et al. 2018), could be used as a preliminary tool to detect outlier specimens. If the sample size of fishes and number of species is sufficiently large, then it could be reasonable to assume that within- and between-species variability is well represented. In this case, LDA (Wakefield et al. 2014) is a useful tool to identify specimens that do not belong to the species list in the analysis with an estimated probability.

The authors declare that they have no conflicts of interest.

This study was based on results of two projects, Fisheries Program I: Population Dynamics and Stock Assessment (Project ID: CEW02374) and Fisheries Program II: Essential Fish Habitats (Project ID: CEW02375), which were funded by the Saudi Arabian Oil Company and run by Marine Studies Section, King Fahd University of Petroleum and Minerals.