SC/F14/J29

Dynamic population segregation bymorphometrics in Antarctic minke whales

genetics

andToshihide Kitakado1, Tore Schweder2, Naohisa Kanda3, Luis A. Pastene3 and Lars Walløe2
1Tokyo University of Marine Science and Technology, 4-5-7 Konan, Minato-ku, Tokyo 108-0075, Japan
University of Oslo, P.O. Box 1103, 0317 Blinder, Oslo, Norway
3
Institute of Cetacean Research, 4-5 Toyomi-cho, Chuo-ku, Tokyo, 104-0055, Japan
Contact e-mail: kitakado@kaiyodai.ac.jp
2

ABSTRACTAn integrated approach for estimating longitudinal segregation of two populations using different sourcesof information, genetic and morphometric data is introduced. A soft boundary is allowed to vary by year
and sex although baseline populations are assumed. A joint likelihood function derived from the twosources is defined for the estimation of mixing proportions. The method is applied to the extensive data
for Antarctic minke whales taken by the JARPA and JARPA II surveys during the austral summers from
1989/90 to 2010/2011 in Antarctic areas III-E, IV, V and VI-W. The mixing proportion is modeled by alinear logistic model with population-specific parameters for the two sets of data. The genetic andmorphometric data in the western areas (Areas III-E and IV-W) and the eastern areas (Areas V-E andVI-W) are used only for estimating parameters for the two putative populations, I-stock and P-stock,
respectively. It was observed that the morphometric data shown dominated information compared to thegenetic data and it helped convergence in the optimization. However, the inclusion of the morphometricdata altered the estimation results and tended to give softer boundaries. On the whole, the result indicatesthat the spatial distribution of the two populations has a soft boundary in Area IV-E and V-W, which
depends on the year. It also suggested possible sex differences along the boundary.

INTRODUCTION
Minke whales are feeding in large numbers in the circumpolar Southern sea during the austral summer. Recently,
the abundance levels of minke whales in the Antarctic Ocean were estimated based on the data derived from the
IDCR-SOWER program (Okamura and Kitakado 2012, Bravington and Hedley 2012, and IWC 2013). In
addition, the estimates were made for areas III-E, IV, V and VI-W from line transect surveys carried out in the
JARPA program (Hakamada and Matsuoka, 2014). These surveys show larger variation in abundance estimates
than what is consistent with the nominal standard deviations representing sampling variability. This might
partially be due to shifting oceanographic and/or feeding conditions in the areas from year to year. This mightalso have an impact on the yearly changes in the underlying stock structure.
The population structure for Antarctic minke whales has been investigated by using genetic and morphometricdata obtained by the JARPA sampling surveys from 1987/88 to 2004/2005. Pastene et al. (2006) found genetic
separation by longitude although the genetic divergence between populations is small. The paper also showedthat Area V-W is a mixing area of two stocks but the transition area may change by year. In addition, a cluster
analysis on a subset of the body measurements also indicates that there are two populations that feed in thisregion of the Antarctic Ocean (Hakamada 2006). These two papers contributed to a concluding paper, Pastene(2006).
Schweder et al. (2011) and Kitakado et al. (2012) also investigated the population structure through an integrated1method by simultaneously using the two sources of information on genetic and morphometric data for the
JARPA data and these papers indicated that there is a soft boundary in Areas IV-E and V-W between the twodifferent populations and the boundary significantly depends on the year. The result also suggested that theboundary might be sex-specific although it showed greater variations from year to year.
Now additional extensive data on genetic and body measurements have also been collected from harvested
Antarctic minke whales from the JARPA II sampling surveys. These data would not only expand the range ofyears to monitor the yearly changes but also contribute to obtaining information on baseline populations for the
mixture analyses.
In this paper, we extend our previous analyses by updating the data sets from the JARPA and JARPA II. Genetic
information was separately used, but the results with morphometric data were also provided. We assume twobreeding populations softly separated by longitude and other covariates according to a linear logistic mixing
model with year-specific intercept.
MATERIALS
Microsatellite data were obtained from analyzing up to twelve sets of primers EV1, EV104, GT023, GT211,
GT195, DlrFCB14, AC045, AC082, AC087, AC137, CA234 and GT129. More detailed information on thelaboratory procedures are described in Kanda et al. (2014). Also up to 10 measurements of body length and ofother lengths between points on the body of the minke whale were available. Each sampled individual isrecorded by date and location of sampling as well as biological information on sex and maturity. Figure 1 shows
sampling locations of minke whales in 1989/90 - 2010/2011 seasons. Table 1 summarizes the sample size by year.
The samples obtained in 1987/88 and 88/89 were not used because of lack of microsatellites.
Here, we denote genetic markers on up to 12 microsatellite loci and the morphometric measurements denote as{(a2l 1 , a2l ); l  1, 2,....,12} and v   v1 , , v10  . Measurements of body length and of other lengths between
points on the body of the minke whale, as shown on Figure 2, were transformed into the 9 allometric measuresagainst the body length as mi  log(vi 1 / v1 ) (i  1, 2,...,9) . The plots of the morphometric data against thelongitude for mature individuals shown in Figure 3 suggest that Eastern individuals tend to have higher values inall but a couple of the measurements. To take into account the allometric change with growth, observed maturity
status for each individal is used.
STATISTICAL MODELLING
Population mixture models
Consider the case of at most two different populations of Antarctic minke whales. Let y represent year andx longitude of a sampled individual. Longitude is here measured in degrees with origo at 180 degrees (E or W),
and with negative numbers to the west and positive to the east. With x  0 the whale was taken at180  x degrees east, and for a positive x it was taken at 180  x degrees west. For example, x= -20 and 20
mean 160E and 160W, respectively.
The probability that an individual belongs to the Eastern population (P-stock in Pastene 2006) is assumed tofollow the linear logistic form,PI  E  p   xe y.
  x1 e y(1)Here   0, and M 50 y   y /  is defined as the longitudinal point of 50% mixing by year. Of course, theprobability that an individual belongs to the Western population (I-stock) is P  I  W   1  P  I  E  . The
2

model with different mixing proportions by year and sex is straightforward. Table 2 summarizes a list ofrepresentative models employed in this study.
To see the effectiveness of such longitudinal segregation, we also estimate area-wise mixing proportions by year
without assuming the longitudinal covariate.
Likelihood components: overview
Each individual contributes likelihoods from its microsatellite and its body measurements, when available. The
parameters in the likelihood components depend on population indexed by i=E,W.
The microsatellite alleles at locus l , G2l 1 , G2l , are assumed independent and identically multinomiallydistributed with probabilities P  G  a    lai for the set of alleles a observed at the locus, and independent
across loci. The 9 morphometric measurements M are assumed multivariate normally distributed with meaniivector  depending on both the population and sex and covariance matrix  depending on sex or
population and sex. The variables G and M are further assumed independent within an individual, and also
between individuals. Note that the morphometric data from immature animals are not used for the analysis, andtherefore the estimated parameters for the immature animals based on both data sets are usually the same as
those from the genetic data unless parameters common to the groups (mature males, mature females, andimmature animals) are incorporated.
Likelihood contribution from microsatellites
Disregarding genotyping errors, and dropping the index for population, a homozygote (a, a ) in a locus l has
2probability gl ,aa   l ,a , while a heterozygote (a, b) has probability gl ,ab  2 l ,a l ,b . The likelihood of theobserved microsatellite alleles of an individual is then
12

12

l 1

l 1

Lmicro  1  p   glW,a2 l1a2 l  p glE,a2 l1a2 l .

(2)Likelihood contribution from Morphometrics
Within population i and sex s, the morphometric measurements are assumed multi-normally distributed withmean depending on population and sex and with covariance matrix only depending on sex or possibly also on
population. With the covariance matrix not depending on population, the likelihood of the vector m of the 9observed measures of an individual whale of sex s is


 1

1/2
1
1  p  |  s | exp   2  m   sW  '  s  m  sW   
9/2

 .
Lmorph   2  


 1

1/2
1
  p |  s | exp    m   sE  '  s  m   sE  

 2




(3)

Joint likelihood
The likelihood of an observation for each individual, conditional on population, is by independence

Li  Lmicro ( i ,  ) Lmorph (Male,i , Female,i , i ),i  E, W .

(4)

The mixing over populations is done at the individual level. Conditioning on “year”, “longitude” and “group”(mature male, mature female and immature animals), and other covariates the likelihood for the observed3individual isL  p( y ,  ) LE  1  p( y ,  ) LW .

(5)

Here, we assume that western areas (Areas III-E and IV-W) and eastern areas (Areas V-E and VI-W) arerespectively occupied by baseline populations (I- and P-stocks). Therefore, the genetic and morphometric data inthese two areas are used only for estimating parameters for the two putative populations. A reason of employing
this baseline assumption is to provide stabilized baseline populations themselves to reduce the possibility ofoccasional disturbance by third or even fourth populations visiting our study areas. A sensitivity test was
conducted to examine the influence of the baseline assumption by assuming narrower baseline areas (Area III-E
and VI-W for the two populations, respectively).
The method is then to maximize the product of all the individual likelihood components. This scheme is carriedout by using the software ADMB for various models for mixing.

RESULTS
Table 2 summarizes the information on model specifications and model selection. The number of parametersmeans the effective mixing parameters estimable by the data. The result indicates that the inclusion of
longitudinal covariate tends to improve the model fitting.
In the estimation based only on the genetic data, the best model is Model (M50$1 + Slope$1), which does notpossess any “year” covariate. Neither the “year” nor “group” effect was statistically significant. In fact, Figure 6ashowed the estimated longitudinal mixing under Model (M50$Year + Slope$1), which indicated less variable
mixing patterns across years. However, this might be attributed to a non-significant masked effect by “group”.
Figure 6b showed “group”-specific yearly variation, which indicated a sharper mixing slope than the combined
one and potentially different mixing patterns among groups.
The “year” effect became clearer when the morphometric data were employed (see Table 2). The result indicatedthat the morphometric data showed dominated information compared to the genetic data and it helped
convergence in the optimization. Meanwhile, as shown in Figures 7a and 7b, the inclusion of the morphometricdata more or less altered the estimation results and tended to give softer boundaries.
Figures 8 and 9 showed the pattern of yearly changes in the longitude at 50% mixture (M50). The parameter
M50 ranged over Areas IV-E to V-W centering around x = 130E to 140E (-50 to -40 in our definition).
Model selection by AIC does not clearly say, but on the whole, the result indicates that the spatial distribution of
the two populations has a soft boundary in Area IV-E and V-W, which depends on the year. It also suggested
possible sex differences along the boundary.

DISCUSSION
We have introduced our integrated analyses for estimating longitudinal segregation of two populations using
genetic and morphometric data. It was demonstrated that the spatial distribution of the two populations has a soft
boundary between E and W populations (so-called I- and P-stocks) and it may depend on the year. The results
also suggested that the boundary may be sex-specific. Some further covariates including environmental
indicators might improve the fitness to the data.
The reason that we integrated these sorts of information for inferring population structure is to use all theavailable data relevant to the structure. More correctly, we originally thought the morphometric data would workby supplementing the genetic data. However, as shown in our previous paper, for the case of Antarctic minke
whales, the likelihood contribution from the morphometric data is much greater than that from the genetic data.
This was the case in this updated analysis. This is attributed to the populations being weakly differentiated. The

4

situation of weakly differentiation is very common for fishery populations, and therefore, we believe ourmethods are of course promising to investigate cetacean population structures but also useful for other species.
As shown in Figures 6 and 7, the estimated smooth mixture functions may not be compatible with our baseline
assumption (Areas III-E and IV-W for I-stock and Areas V-E and VI-W for P-stock). To see the impact of thebaseline assumption on the estimation result, a sensitivity test with a different baseline assumption (Areas III-E
and VI-W are baseline line areas for I- and P-stock, respectively) was conducted with both the genetic andmorphometric data (see Figure 10). Some difference in the estimation results was observed, which means that aslight sensitivity to the baseline assumption may occur. It is hard to confirm that the assumed baseline is the truebaseline or not, but more reliable analysis with some more valid baseline assumptions may become possible if
the samples outside the JARPA and JARPAII areas are available. This warrants further investigation with abroader longitudinal range of data.
For handling morphometric data of immature animals, a model expressing allometric change along with growth
can be developed (see Figure 11). This model has three parameters; a mean of allometric quantity in the maturestate, a logarithm of length at which the shift of phases of allometic relationship between immature and mature
occurs, and a slope parameter expressing allometric changes in the immature state. We noted that the second
parameter relevant to the shift in phase is shared with all the (nine) measurements, which means the parameter isassumed to be common to the measurements: a  b (log v1  c)
E[mi ]   i i
ai

for log v1  cfor log v1  c(6)

The model considers a change in allometric relation between immature and mature stages. The allometric modelis originally derived by the famous formulas;vi   Imm v1Imm

for immature animalvi   Mat v1Mat

for mature animal(7)

with the assumption of Mat  1 , for which no allometric change is supposed to happen in the mature state, andthe connectivity at the change point. The parameter  Mat might not necessarily be 1, so we will also be able
investigate the following general case:ai  bImm,i (log v1  c)
E[mi ]   ai  bMat,i (log v1  c)for immature animal
for mature animal(8)

Furthermore, the change point is not necessarily common to the measurements but may differ between spine andcranial bones, and hence this is also to be tested in the model. Unfortunately, this investigation has not been
completed yet because of missed communication in the data set between the data owner and analysts, but thismodelling is highly worth conducting.
As mentioned earlier, microsatellites could not be read for some individuals at some loci. We investigated several
models for such drop out. The question of whether this is caused by allelic drop out is of specific interested. The
presence of allelic drop out would weaken the information content in the recorded microsatellite data, and thus
affect its likelihood function. If also the probability of allelic drop out depends on population, the estimatedprobabilities of population assignment would be biased. The size of such biases could potentially be investigatedby simulation, something which is not done here. Under individual and locus drop out, i.e. all loci dropping out,
Hardy-Weinberg equilibrium can be tested by contrasting the number of homozygotes to its expected value under
the null hypothesis of equilibrium. Consider one locus l . The number of homozygotes recorded for this locus
is Tl . The estimated probability for a given individual to be homozygous is the mixture of the estimatedE
population specific probabilities of homozygotisity for an individual i , qi ( pˆ i gˆ aa
 (1  pˆ i ) gˆ Waa ) .
Individuals are independent with respect to homozygosity, and E (Tl ) qi anda var(Tl )  qi (1  qi ) .
There is also independence between loci. The mean and variance under thei null hypothesis of thei test statistic



5







T  Tl are thus obtained by summing the locus specific quantities. The normal approximation should be
excellent in our large data set, and Hardy-Weinberg equilibrium is tested by homozygote excess in the usual way.
This test might be regarded as a conditional test since the data for estimating the null distribution are ancillary
to T . This kind of exercise is also under way.

REFERENCES
Bravington, M.V. and Hedley, S.L. (2012) Abundance estimates of Antarctic minke whales from the IWC
IDCR/SOWER surveys, 1986-2002. Paper SC/64/IA13
Hakamada, T. (2006) Morphometric analysis on stock structure in the Antarctic minke whale based on JARPA
samples. SC/D06/J10
Hakamada, T. and Matsuoka, K. (2014) Estimates of abundance and abundance trend of the Antarctic minke
whale in Areas IIIE-VIW, south of 60S, based on JARPA and JARPAII sighting data (1989/90-2008/09).
SC/F14/J3
International Whaling Commission. (2013) Report of the Scientific Committee. Annex G. Report of the
Sub-Committee on In-Depth Assessments. J. Cetacean Res. Manage. (Suppl.) 14:195-213.
Kanda, N, Goto, M, Oikawa, H. and Pastene, L.A. (2014) A note on sampling and laboratory procedure protocolsof the genetic work at the Institute of Cetacean Research. SC/F14/J27
Kitakado, T., Schweder, T., Kanda, N., Pastene, L.A. and Walløe, L. (2012) Progress report on the estimation of
longitudinal mixing proportions for the Antarctic minke whales using genetic and morphometricmeasurements. SC/64/IA4
Okamura, H. and Kitakado, T. (2012) Abundance estimates of Antarctic minke whales using the OK method.
Paper SC/64/IA2
Pastene. LA. (2006) What do we know about the stock structure of the Antarctic minke whale? A summary ofstudies and hypotheses. SC/D06/J12.
Pastene, LA., Goto, M. and Kanda, N. (2006) Genetic analysis on stock structure in the Antarctic minke whales
from the JARPA research area based on mitochondrial DNA and microsatellites, SC/D06/J9.
Schweder, T., Kitakado, T., Kanda, N., Pastene, L.A. and Walløe, L. (2011) Dynamic population segregation by
genetics and morphometrics in Antarctic minke whales. SC/63/IA7.

6

Table 1. Sample size by year and area.
Year
1989/90
1990/91
1991/92
1992/93
1993/94
1994/95
1995/96
1996/97
1997/98
1998/99
1999/00
2000/01
2001/02
2002/03
2003/04
2004/05
2005/06
2006/07
2007/08
2008/09
2009/102010/11

Mature male
IV-W IV-E V-W V-E VI-W
114
42
99
48
96
32
105
45
72
77
97
60
47
130
30
38
70
52
60
88
52
83
92
26
39
47
81
71
71
49
43
53
50
48
83
45
53
59
33
33
58
52
48
216
27
73
74
51
107
54
6
28
55
81
157
81
33
76
2
16
25

Total

478

III-E

962

430

806

700

457

Mature female
IV-W IV-E V-W V-E VI-W
45
40
57
74
61
23
66
57
38
31
14
73
29
65
11
58
113
22
19
23
5
62
18
11
28
73
26
85
17
31
43
76
36
84
11
27
106
36
43
143
18
37
153
13
38
1
272
16
32
95
3
52
31
130
39
113
49
0
37
10
84

III-E

299

706

311

7

520

1060

207

Immature animals
IV-W IV-E V-W V-E VI-W
63
22
31
14
53
23
34
20
47
65
29
57
33
62
33
31
34
22
32
85
74
38
35
35
59
25
76
43
34
44
36
58
50
16
64
53
30
45
51
6
47
40
45
130
34
37
1
58
34
90
73
4
7
48
39
99
51
21
35
8
20
15

III-E

376

662

432

355

431

342

Table 2. Description and comparison of models employed in this study. The highlighted models are not the bestmodel in terms of AIC, but those are used for the graphical presentation of results in Figures 6 to 9. The numberof parameters is that only for estimable parameters in the mixing proportions. Convergence was not reached insome cases due to an inadequate gradient and/or degenerated hessian matrix.
Data
Genetic

Genetic

Both

Both

Convergence

negloglike

Delta-negloglike

#para

AIC

Delta-AIC

Area
1
Group
Year
Year*Group

OK
OK
OK
Grad OK. Hess No

294734294733294724294711

40
39
30
17

2
6
21
66

589472589478589490589554

6
12
24
88

M50
Slope
1
1
Group
1
Group
Group
Year
1
Year
Group
Year
Year
Year*Group1
Year*GroupSM
Year*Group Year*Group

OK
OK
OK
OK
OK
Hess OK, Grad No
NO
NO
NO

294731294731294730294720294715294713294694294715294732

37
37
36
26
21
19
0
21
38

2
4
6
21
25
44
67
69
132

589466589470589472589482589480589514589522589568589728

0
4
6
16
14
48
56
102
262

OK
OK
OK
OK

98465.8
98463
98439.898416.3

55.7
52.9
29.7
6.2

2
6
21
66

115.4
117.8
101.4
144.4

26
28.4
12
55

OK
OK
OK
OK
OK
Hess OK, Grad No
OK
No
No

98455.3
98453.998452.898432.1
98429.8
98428.7
98410.1
98430.2
98449.8

45.2
43.8
42.7
22
19.7
18.6
0
20.1
39.7

2
4
6
23
25
44
67
69
132

94.4
95.6
97.4
90
89.4
125.2
134
178.2
343.4

5
6.2
8
0.6
0
35.8
44.6
88.8
254

Area
1
Group
Year
Year*Group
M50
Slope
1
1
Group
1
Group
Group
Year
1
Year
Group
Year
Year
Year*Group1
Year*GroupSM
Year*Group Year*Group8

Figure 1. Locations of sampling by year. Blue, red and green dots show mature males, mature females andimmature animals. The values on the horizontal axis showed the longitude (e.g. -100, -50 and 0 respectively
mean 80E, 130E and 180E).

9SC/F14/J29
v1: Body lengthv2: from the tip of snout to center of eye
v3: from the tip of snout to ear
v4: from the tip of snout to tip of flipper
v5: from the tip of snout to end of ventral gloves
v6: from the tip of snout to center of umbilicus
v7: from the tip of snout to sexual apparatus
v8: from the tip of snout to anus
v9: length of skull
v10: width of skull

Figure 2. Morphometric measurements for Antarctic minke whales used in this study. These measurements otherthan v1 are transformed to the logarithms of allometric measures asmi  log(vi 1 / v1 ) .

Figure 3. Morphometric data against the longitude where samples were taken (The blue and red dots show mature
male and mature female individulas, respectively).

10(a) Based only on the genetic data(b) Based on the genetic and morphometric data

Figure 4. Yearly variation of mixing proportion by Area and Group (blue=mature males, red=mature females,
green=immature animals, black=all animals).

11Figure 5. Comparison of area-wise mixing proportions estimated by the genetic data and by the genetic andmorphometric data (genetic data only=solid line with open circles; genetic and morphometric data=shaded line withfilled circles) for mature males and females as well as immature animals and all individuals ((blue=mature males,
red=mature females, green=immature animals, black=all animals).

12Figure 6a. Estimated mixing proportions against longitude for the Eastern population (P-stock) based only on thegenetic data under Model (M50$Year + Slope$1) that the intercepts in the mixing proportions differ across years
while the slopes are common to years and groups. Circles are estimated mixing proportions where samples weretaken.

13

Figure 6b. Estimated mixing proportions against longitude for the Eastern population (P-stock) based only on thegenetic data under Model (M50$Year*Group + Slope$1) that the intercepts in the mixing proportions differ acrossyears and groups while the slopes are common to years and groups. The blue, red, and green lines are for mature
male, mature female and immature animals, respectively. Circles are estimated mixing proportions where samples
were taken.

14Figure 7a. Estimated mixing proportions against longitude for the Eastern population (P-stock) based on the genetic
and morphometric data under Model (M50$Year + Slope$1) that the intercepts in the mixing proportions differ
across years while the slopes are common to years and groups. Circles are estimated mixing proportions where
samples were taken.

15Figure 7b. Estimated mixing proportions against longitude for the Eastern population (P-stock) based on both thegenetic and morphometric data under Model (M50$Year*Group + Slope$1) that the intercepts in the mixing
proportions differ across years and groups while the slopes are common to years and groups. The blue, red, andgreen lines are for mature male, mature female and immature animals, respectively. Circles are estimated mixing
proportions where samples were taken.

16Figure 8. Yearly variation of longitudes at which 50% mixtures occur. The assumed model is Model(M50$Year*Group + Slope$1). The blue, red, and green dots and lines are for mature male, mature female andimmature animals, respectively. The black ones are the estimates from Model (M50$Year + Slope$1).

17Figure 9. Yearly variation of longitudes at which 50% mixtures occur by group. The assumed models are same as in
the Figure 8.

18

Figure 10. Estimated mixing proportions against longitude for the Eastern population (P-stock) based on both thegenetic and morphometric data under Model (M50$Year*Group + Slope$1) with a different baseline assumption(III-E and VI-W are respective baseline populations).

19

mi  log(vi 1 / v1 )ai

ci

log v1Figure11. An assumed allometric change in the immature period. (c’s are common to measurements)

20