Approaches to otolith age determination : image signal treatment and age attribution *

The dynamic population models used in the management of fishery resources and the biological studies of fish, require age data for the determination of the composition by age for the catch and of the growth rates (Ricker, 1973). The age determination is traditionally carried out from the seasonal growth marks that are formed in the calcified tissues (scales, bones, fin rays, otoliths) (Williams and Bedford, 1974). Given the correlation between the size of the fish and that of the bony structure, it is possible to determine the growth rates, by the measurement of the distances between successive growth structures. The seasonal variations of growth can be studied as a function of the environmental conditions and of the involved genetic stock (Carlander, 1987; Campana, 1987). Since the 1970’s and in parallel to the development of relatively cheap computer equipment, the implementation of semi-automatic measurement programmes and age determination based on the otoliths (Mason, 1974) or the circulii growth rings in scales were initiated (Troadec and Prouzet, 1986; Gandelin and Laval, 1987). The discovery that growth structures were formed with a daily periodicity in the otoliths (Pannella, 1971), exponentially increased the use of the daily growth increments in the study of the growth variations (Gutiérrez and Morales-Nin, 1986), recruitment (Victor, 1982) and the relationship environmental parameters and survival (Stevenson and Campana, 1992). Evidently, the problems of age SCI. MAR., 62 (3): 247-256 SCIENTIA MARINA 1998


INTRODUCTION
The dynamic population models used in the management of fishery resources and the biological studies of fish, require age data for the determination of the composition by age for the catch and of the growth rates (Ricker, 1973).The age determination is traditionally carried out from the seasonal growth marks that are formed in the calcified tissues (scales, bones, fin rays, otoliths) (Williams and Bedford, 1974).Given the correlation between the size of the fish and that of the bony structure, it is possible to determine the growth rates, by the measurement of the distances between successive growth structures.The seasonal variations of growth can be studied as a function of the environmental conditions and of the involved genetic stock (Carlander, 1987;Campana, 1987).
Since the 1970's and in parallel to the development of relatively cheap computer equipment, the implementation of semi-automatic measurement programmes and age determination based on the otoliths (Mason, 1974) or the circulii growth rings in scales were initiated (Troadec and Prouzet, 1986;Gandelin and Laval, 1987).
The discovery that growth structures were formed with a daily periodicity in the otoliths (Pannella, 1971), exponentially increased the use of the daily growth increments in the study of the growth variations (Gutiérrez and Morales-Nin, 1986), recruitment (Victor, 1982) and the relationship environmental parameters and survival (Stevenson and Campana, 1992).Evidently, the problems of age SCI. MAR., determination increased when growth structures with an annual periodicity went to daily structures.In various laboratories semi-automatic methods based on the use of a TV camera coupled to a microscope and a computer were implemented (Ralston and Williams, 1989).The complex growth of the otolith and errors in its preparation for observation, frequently make the resolution of the growth structures difficult.On the other hand in certain phases of the growth the formation of rings with a periodicity of less than 24 h is frequent.Image analysis is a new tool that allows the images to be improved and makes interpretation easier (Campana, 1992).
Another characteristic of the otoliths is that their shape and relative size are specific to each species (Schmitt, 1969;Gaemers, 1976;Nolf, 1985, Lombarte andMorales-Nin, 1995), implying a genetic regulation in their formation.Consequently, the morphological study of the otoliths can be very useful in the separation of species and in the determination of phylogenetic lines (Gaemers, 1984;Nolf and Steurbaut, 1989;Lombarte and Castellón, 1991).The ultrastructural and morphological study of the otolith has also brought knowledge about their function in fishes inner ear (Gauldie, 1988;Lombarte, 1992).
As there is a relationship between "gross" otolith morphology and its ultrastructure (Gauldie and Nelson, 1990;Lombarte and Morales-Nin, 1995), a tool for the otolith study should contemplate both aspects.
However, the development of algorithms for image analysis applied to ageing is still ongoing (EFAN, Rep. 1997).In the present study some new approaches in age applied to an image analysis system are discussed in an attempt to provide new tools for this methodology.Although the method and discussion refers to increments, it can be applied to annual zones.

Image obtention and signal treatment
The PC based image systems are used in order to determine the number of growth structures and the relative distances between them.Once the otolith is identified, the procedure to be followed in order to read it is to digitise the image and measure the otolith radius selected for the age determination.
Once the reading radius is measured, the rings that are found throughout it should be detected.For this the image of the otolith is digitised at a magnification at which the increments are visible (Fig. 1).Given that, as a rule, this magnification never shows the complete otolith, the detection of the rings of the reading radius should be carried out by stages, digitising in each one of them the visible zone of the otolith and reading the increments that are found in the segment corresponding to the reading radius.
From the digitised image the profile of the segment is obtained (Fig. 2).The profile of a vector of any digital image is an integer value vector which corresponds with the grey levels of the points of the image that fall below this vector.These grey levels are integer numbers that represent the brightness of the image at a determined point.Its values vary between 0 (level of minimum brightness, equivalent to black) and 255 (level of maximum brightness, equivalent to white).
We start, then, with a N-length vector or signal, x(t).As the rings of the otolith are dark bands, in principle, the problem comes down to finding the local minimum (valleys) of the signal, that is to say, points of the signal whose grey level is less than those to the left and right.
However, as we can observe in Fig. 2a, not all the signal minimum correspond with increments of the otolith.Some of these minimums are a consequence of the noise introduced in the signal in the process of image digitisation.Others, can be sub-rings formed to a greater frequency than that of the studied rings.On the other hand, the radial structure of the aragonite microcrystals can be superimposed on the signal, introducing error.This means that we have to "treat" the signal at various levels before being able to detect the minimum that are the object of our interest: 1) At the level of image capture.In order to minimise the noise coming from the digitisation process, it is recommended that the definitive digitised image is the mean of a series of digitisation's of the same analogical image.Once the image is digitised, we can make use of the utilities of the image analysis system in order to improve its quality: enhance the contrast, optimise the overall illumination, etc.
2) At the level of obtaining the profile.Instead of obtaining the signal uniquely and exclusively from the reading radius profile, it should be obtained as the mean from the profiles of k adjacent segments, parallel to and centred on the reading radius; that is to say, as if the profile were obtained from a radius of k pixels in thickness.These two steps are done before the processing of the signal.From this moment knowledge of the theory of analysis and processing of periodic signals can be applied.

Fourier analysis
The periodic signals can be represented in two ways, basically: 1) Representation in the temporal domain: x(t).In our case, each instant t= 0, 1, ... corresponds with the point t of the segment of the reading radius and x(t) is the grey level (0...255) of this point.
2) Representation in the frequencies domain: The basic concept of this representation is that it considers the signal to be the sum of a series of sinusoidal signals of distinct periodicity, amplitude and phase.In this way, the representation of the signal in the frequencies domain, X(u), shows us for each frequency u, the amplitude and the phase of the corresponding sinusoidal signal.
Given the nature of the problem, we consider that the optimum representation of the signal is the second.The method used in order to obtain the representation of the frequencies in space for the signal x(t) is the Digital Fourier Transform (DFT), implemented in an algorithm called FFT (Fast Fourier Transform).The result is two real valued vectors xreal (u) and ximag(u), that correspond to the real part and the imaginary part, respectively, for a complex number.If we represent this complex number in its polar form (model, angle), the amplitude of the sinusoidal signal of frequency u, it would correspond with the module and the phase with the angle.
In the ideal case that the signal x(t) would be composed by a simple sinusoidal signal of frequency k and amplitude A k , the number of rings (dark bands) would be, precisely, k.But, as we have previously seen, we found signals that, by themselves, are not composed of an unique signal, but of multiple signals of various frequencies.
As a general rule, the signals of higher frequencies correspond to the original signal noise, coming from the image digitization process.Consequently, in order to eliminate this noise, we filter the high frequency signals the original signal is composed of, onwards from a moving threshold.That is to say, we cancelled out the value of all the signals of frequencies higher to one given.This threshold cannot be fixed, since, depending on the species and the growth phase, the otolith can show increments of various thickness.Another factor that influences the number of rings by image is the magnification used.For all this, we can find both an image in which the rings are very fine and are very close to one another (high frequency), and one in which the rings are thicker and are more separated (low frequency).Therefore, this threshold has to be established for each reading.
Once the frequencies corresponding to the noise are eliminated, in order to be able to detect the exact position of each ring, we should have its representation in the temporal space.For this the inverse of the Digital Fourier Transform is calculated (that leads to the algorithm FFT being applied again to the combined X´(u)), that is to say, to the vectors xreal´(u) and -ximag´(u), so obtaining a filtered version of the signal x(t) : x´(t).
The local minimums of this filtered signal, x´(t), will indicate to us the position of the rings.The signal that we can observe in Fig. 2b has been obtained from the signal of Fig. 2a, resulting from eliminating all the frequencies greater than 13.If we find the local minimums of this new signal, we locate the 10 increments that are observed in the reading radius.

Wavelet analysis
Generally the signals are studied in the frequential domain instead of the spatial domain means of Fourier analysis.The transformed Fourier Function measures the high frequencies of the signal but does not locate their position.This problem has been solved by the wavelet analysis developed from the functions ψ (x), which transformed, ψ (sx+t) can be used to expand functions .This results in an intermediate representation between the spatial and frequential domains.Thus, the irregularities of the function can be placed and the structures found in the signal at different levels can be studied with an increasing resolution.
Thus, the signal function ψ (x) is decomposed in a family of functions resulting from the translation and dilatation of the function.This function ψ is the ondelette, the corresponding family is: The transformed of the ondelette function is: If the dilatation of ψ is the factor s, then: Thus the transformed ondelette is the following scalar product: which corresponds to the decomposition of f in a family of functions: ψ s (x) has the same shape as ψ(x), but with a size S times smaller (Fig. 4).

If:
the transformed ondelettes in a point u with scale s can be: Wf (s,u) = f(x) ψ - s (u) this is the filter of f by the band-pass filter ψ - s (x).These transformations can be applied to the signal from the reading vector analysing together the temporal and frequency domains.
The resolution r of the signal gives the minimum size of the structures that can be detected by the signal.

Age attribution
During the reading a weight is given to every ring, in the following manner: -if the ring has been detected by the algorithm or by the user, it's given a weight of 1.
-If the ring has been obtained as a result of an interpolation, the weight will be given by the ratio detected rings/(detected rings + interpolated rings).
In order to minimize possible errors, it is advised that each otolith is read at least twice.For each of the TR complete readings of an otolith (k:1..TR) we obtain the following data: the measurement of the reading radius (R k ), the number of rings (X K ) and, for each ring (j), its distance to the centre of the nucleus (r kj ) and a factor of quality (w kj ).
Once the TR readings have been carried out on the same otolith, the data of the various readings should be combined into one.
In order to be able to compare the data from the various readings, these should be standardised.The larger radius among all the readings is selected (R max ) and, for each one of them, the distances of the read rings are multiplied by a scale factor: larger radius/reading radius.
In order to obtain the age estimation we follow the procedure proposed by Campana (1992) with some changes.Campana supports that the best age estimation is obtained from multiple readings of each of several otoliths of the same fish.An appropiate age estimation procedure would involve: 1-the determination of the single best estimate of increment count for each otolith: ; Where k is the reading index t is the otolith type (sagitta, lapillus or asteriscus) j is the otolith (left or right) TR tj is the number of times the otolith has been read W tjk is the otolith reading weight X tjk is the increment count for a reading The weight W tjk is a value in an arbitrary scale, from 1 (little confidence) to 5 (unambiguous count), on the basis of perceived confidence in the reading.
2-pooling of the increment counts from a given otolith type.

;
Where t is the otolith type (sagitta, lapillus or asteriscus) r is the right otolith l is the left otolith 3-converting increment counts from each otolith type to age estimates.
The necessary data to realise this conversion are: age of the fish at the first increment deposition (A t0 ) and the increment formation frequency, f t − ψ (days/increment).As a general rule, we'll consider daily increments: f t =1.
A t = A t0 +f t C t 4-pooling the results from the otolith types.Our approach differs from the one above only in the first step, when determining the best count for each otolith: instead of calculating the median of counts, we'll try to fuse rings from different readings that refer to the same otolith increment.For instance, if we have the following results from three readings of the same otolith (see Table 1).Following Campana's procedure (1992), the best otolith count would be seven rings, because in the three readings seven rings have been counted; but if we represent graphically these readings (Fig. 3), we can easily see that it would be more convenient to find out which rings from different readings actually refer to the same increment in the otolith, and obtain its position as a combination of the positions of the related rings.But, how can we determine which rings from different readings refer to the same otolith increment, that is to say, are combinable?We consider the following heuristic rules: -two rings refer to the same increment and consequently are combinable if the distance between them is less than the minimum distance between rings from the same reading (Dmin) (Fig. 4a  -if a ring can be combined in several ways the best combination will be the one that combines closest rings (Fig. 4b).
As a result of combining the three readings of the example above we obtain a new sequence of eight rings (see Table 2).
The algorithm for a combination of rings consists of a series of iterations in which each one tries to combine XC rings.At the beginning, XC is equal to TR, that is to say, to the total number of readings.When this number of rings cannot be combined any more, the next step is to decrease XC by 1 unit, and try the combination again in this way, until all possible pairs of rings have been combined.
For each increment we have a set of selected rings from different readings that are to be combined: In order to find the distance to the centre of the nucleus from the ring resulting from the combination of several readings, a weighted mean is calculated from the various distances.The weighting is given by the quality factor of the reading.

;
Where: cr i is the distance from the resulting ring to the centre of the nucleus cw i is the weighted quality of the combination of the ring sr im is the distance to the centre of the nucleus from the m th ring to be combined sw im is the corresponding quality XC is the number of readings that are to be combined TR is the total number of readings Once all the otoliths of the individuals from the same species have been read, the data are summarised in the form of a matrix in which the rows correspond to the studied individuals and the columns to the data of each individual: sample number, weight, size, sex, length in mm for the reading radius, number of increments read, date of capture and date of birth.The date of birth is obtained from the following data: date of capture, number of incre-  ments, frequency (in days) of formation for each increment and the number of days that the individuals of this species take to form the first increment.Generally, the frequency of increment formation will be daily.In case the duration of the period of formation of the first increment is unknown, it is supposed that it is 0 days.

DISCUSSION
Since the pioneering work of Mason (1974) based on seasonal growth structures, numerous methods of age determination based on semiautomatic systems have been developed (Method 1981;Frie, 1982;Messieh and McDougall, 1985;McGowan et al., 1987;Tzeng and Yu, 1988;Ralston and Williams, 1989;Troadec, 1991).They are based on a high resolution video camera coupled to a microscope and connected by means of an interface to a computer.In the monitor of the system the material to study (scale, otolith, bony radius, etc.) is visualised and the rings are read manually, marking both the reading radius and the thickness of the rings, that are registered in an ASCII file.
From this first generation of reading systems, followed a principle of automation based on the light spectrum.The system selects some of the growth structures (Ralston and Miyamoto, 1983).This interpretation is completed by the user in an interactive way.The more recent methods incorporate the image analysis systems in an attempt to improve the quality of the image obtained from the microscope before processing the signal (Campana, 1987;Planes et al., 1991;Chauvelon et al., 1992).In the majority of methods not all the information that is supplied by the digital image (2D) is used, except the information that corresponds to the reading radius (1D : signal or vector).This is due to the fact that the treatment of all the information the image contains would lead to an unacceptable response time.However, the structure of the signal is complex, including biological information (differences in growth patterns due to environment and physiological cues), growth patterns of a different timescale as well as errors due to the preparation.
To reduce the signal complexity and obtain meaningful results a filtering is necessary.Some methods include biological information in order to interpret the growth structures after their automatic acquisition (Small and Hirschhorn, 1987;Planes et al., 1991;Troadec, 1991).For instance, the sup-posed position of the rings can be fitted following a growth model of Von Bertalanffy (Small and Hirschhorn, 1987).Each ring was rejected or accepted taking into account the residuals of the fitting.Also the biological knowledge might consist of a growth pattern, which will take into account both the size variability of the rings and the non linearity of its growth.The number of rings might be estimated by Fourier spectral analysis (the frequency with the greatest amplitude) (Troadec, 1991).The inconvenience of these systems, in which biological knowledge is introduced previously, is their limitation to the study of very specific structures (e.g.otoliths of one certain species of fish), of which the growth pattern is already known.They are not useful for the reading of calcified structures of unknown species or of various growth cases (larval, juvenile, immature, mature).
The signal can be treated directly or better still can be analysed previously by means of temporary series (Bach and Chauvelon, 1992) or other methods (Planes et al., 1991: Troadec, 1991).The vector can be interpreted as a temporary series; the calculated periodogram is used in order to estimate the frequency of the signal and, therefore, the number of periods (micro-increments) that the signal contains.A numerical process of this vector of grey levels is carried out with a moving mean of three points, although other techniques such as spectral analysis and polynomial curve fitting are under development.
We propose the ondelettes, a multiresolution analysis of a signal, as a tool for studying the periodic signals obtained from the otolith reading radius.This method allows the analysis of signals resulting from phenomena produced at different scales.Thus, the analysis considers both the time and frequency domains.This approach is more powerful than simpler statistical treatments of the signal, having the potential to discriminate between the different time-scales (daily, weekly, monthly) interacting with the otolith growth rate.
Our approach synthesises various applications since on the one hand it uses image analysis and on the other hand it uses the luminous signal by means of the Digital Fourier Transformation.We recommend using an interactive approach to change the detection threshold, interpolate with a moving mean, insert nearby zones, or simply mark directly the increments that frequently will be missing due to the lack of a clear growth pattern.On the other hand, the combination of several readings of the same otolith is generally calculated as a mean of the dif-ferent readings.The proposed statistical treatment of the results for several readings of the same otolith, allows an estimate of the considered age to be obtained and sets a criterion for the identification and combination of increments.
Further work on the application of these analytical methods to the vectors of light frequencies obtained from the reading radius is necessary to decod all the information contained in the otolith.Also, it is necessary to consider 2-D and 3-D otolith models in order to take into account all the growth information in the otolith.These developments are feasible thanks to the new powerful and quick computers that open the possibility of short time responses.

FIG. 2
FIG. 2. -a) Signal vector corresponding to a grey level of enhanced increments.b) Top signal with noise and bottom signal once filtered by a Fast Fournier Transform.
FIG.3.-Combination of 3 readings of the same radius (r ij ) in an iterative process (a-c), resulting on a single interpretation (d).

FIG. 4
FIG. 4. -Combination of three readings based on the distances between rings; a, Dmin as limit for combination; b, selection process.

TABLE 1 .
-Example of three different meassurement ridings for a given otolith; r kj radius of the otolith k to the ring j; w kj factor of quality of the otolith k and ring j (see text).

TABLE 2 .
-Results of combinig the reading given on Table1; as in Table1, sr i is the ratius of the ring to bo combined and cw i is the weghted quality of the combination (see text).