A revision of the geographical distributions of the shrews Crocidura tanakae and C. attenuata based on genetic species identification in the mainland of China.

The Taiwanese gray shrew (Crocidura tanakae) and Asian gray shrew (C. attenuata) are so similar in size and morphology that the taxonomic status of the former has changed several times since its description; C. tanakae has also been regarded as an endemic species of Taiwan Island. In recent years, molecular identification has led to several reports of C. tanakae being distributed in the mainland of China. In this study, we determine the geographical distribution of C. attenuata and C. tanakae based on more than one hundred specimens collected during 2000 to 2018 over a wide area covering the traditional ranges of the two species in the mainland of China, and show a substantial revision of their distributions. Among 110 individuals, 33 C. attenuata and 77 C. tanakae were identified by Cytb gene and morphologies. Our results show, (1) C. attenuata and C. tanakae are distributed sympatrically in the mainland of China; (2) contrary to the previous reports, the distribution range of C. attenuata is restricted and much smaller than that of C. tanakae in the mainland of China; (3) Hainan Island, like Taiwan Island, is inhabited by C. tanakae only according to the present data.

Introduction

The Taiwanese gray shrew ( Kuroda, 1938) and Asian gray shrew ( Milne Edwards, 1872) are distinct species with very similar morphological characters and measurements, such that the taxonomic status of has been changed several times by taxonomists. was originally described from Taiwan as a new species by Kuroda (1938); however, because it could not be distinguished from in morphological characters and measurements, was thereafter regarded as a synonym or subspecies, by many authors (Ellerman and Morrison-Scott 1951; Jameson and Jones 1977; Corbet and Hill 1992; Hutterer 1993; Fang et al. 1997;). Motokawa et al. (2001) recognized the distinct taxonomic position of by chromosomal data, and regarded it as the endemic species of Taiwan Island.

In recent years, the application of molecular identification techniques led to reports of populating the mainland of China. Esselstyn et al. (2009) and Esselstyn and Oliveros (2010) genetically identified specimens collected in Vietnam and the Hunan and Guizhou Provinces of China and found most of their specimens belonged to ; only a few were attributed to . Bannikova et al. (2011) and Abramov et al. (2012) reported that was also found in Vietnam and Laos, and it was a wide-spread species in Vietnam, whereas inhabited only the north and east of the Red River; Chinese scientists recently reported was collected from the mainland of China including Mount Emei of Sichuan Province, Mount Fanjing of Guizhou Province, Pingbian and Funing of Yunnan Province and Xingshan of Hubei Province (Cheng et al. 2017; Chen et al. 2018; Lei et al. 2019). However, these reports only provided the data for several distribution areas and were not sufficient to generalise the overall distributions of the two species in the mainland of China. The current IUCN distribution maps of and presented in Figure 1 are revised by this study.

Figure 1.

Distributions of and by the IUCN.

We accumulated more than one hundred specimens from 19 areas of and in our field surveys in the mainland of China from 2000 to 2018, which expands the previous distributions from the aforementioned reports from a few localities. A re-evaluation of geographical distributions of the two species is important to a range of studies and practical needs, such as zoogeography, geophylogeny, agriculture animal management, health and epidemic prevention. Here we report the wide geographical distributions of and in the mainland of China.

Materials and methods

Samples and DNA sequencing

A total of 110 specimens of the complex were collected from 11 provinces in the mainland of China including Anhui (2 specimens), Fujian (4), Guangdong (23), Guangxi (4), Hubei (4), Hunan (3), Jiangxi (13), Sichuan (Chongqing is included) (30), Yunnan (1) and Zhejiang (23) as well as the Hainan Island (2) using Sherman live cages during 18 years from Aug. 2000 to Sep. 2018. The geographical position (longitude and latitude) of each specimen was recorded by GPS (Suppl. material 1, Table S1). All specimens including pelt, carcass and skulls were kept in Shandong University (Weihai) and Guangzhou University.

DNA was extracted from muscle samples using the Ezup Column Animal Genomic DNA Purification Kit (Sangon Biotech Co., Ltd., Shanghai, China). The complete mitochondrial cytochrome b gene (Cytb, 1140 bp) was amplified by PCR with the primers as in Irwin et al. (1991). Primers Nivicob1 (5’-TGTCATTATTTCTACACAGCACTTA-3’) and Nivicob2 (5’-TTTGGGTGTTGATGGTGGG-3’) were used for amplification of the whole Cytb gene. PCR reactions were 25 μL, containing 0.25 μM primers, 2×EasyTap PCR Supermix 12.5 μL, and approximately 15 ng DNA template. The thermocycling protocol was as follows: an initial denaturation of 5 min at 95 °C; 32 cycles of 95 °C for 30 s, annealing temperature (Tm) for 30 s, 72 °C for 1 min; a final extension of 10 min at 72 °C. PCR products were directly sequenced by Sanger sequencing technique.

Phylogenetic analyses

Cytb gene sequences were aligned using BioEdit v.7.2.5 (Hall 1999). Each specimen was molecularly identified for species by blasting on GenBank and confirmed by ML (maximum likelihood) phylogenetic tree construction in MEGA 5 (Tamura et al. 2011) based on TN93+G model. We used the Akaike Information Criterion (AIC) in jModeltest1.0 (Posada 2008) to select the best-fit model of sequence evolution for the locus alignment. The bootstraps were obtained using a rapid bootstrapping algorithm with 1000 replicates. We calculated the genetic distance of Kimura-2-parameter (K2P) of Cytb between the two species.

We also included Cytb sequence data from several earlier studies (Ohdachi et al. 2004, 2006; Bannikova et al. 2006, 2009, 2011; Jenkins et al. 2009, 2013; Esselstyn and Oliveros, 2010; Abramov et al. 2012; Chen et al. 2016) to place the shrews from type locality and Vietnam into a phylogenetic context, the sequence information was showed in Suppl. material 1, Table S2. was selected as outgroup (Suppl. material 1, Table S2). GenBank accession numbers for the original sequences used in this study were MK765682-MK765791 (Suppl. material 1, Table S1).

Morphological analyses

In order to attribute these genetic lineages to taxonomically correct species names, we photographed the dorsal, ventral, lateral of skull and lateral view of the mandible of from type locality – Baoxing (Moupin), Sichuan – and also photographed the corresponding teeth, and marked the characteristic features on the pictures for this species. We repeated the same procedure with the only sample of from the same locality (Baoxing) for interspecific comparisons.

We conducted a morphological investigation of the specimens sampled to identify the two species by determining three external measurements: total body length (TBL), head and body length (HBL), ear length (EL); and 10 skull measurements: greatest length of skull (GLS), cranial base length (GBL), median palatal length (MPL), length of teeth row (LUTR), greatest palatal breadth (GPB), breadth of occipital condyles (BOC), greatest breath of braincase (BBC), interorbital breadth (IOB), height of the braincase (HB), length of mandible (LM) according to Yang et al. (2005, 2007) and Jenkins et al. (2009).The measurements of the skull indices were performed with a digital vernier caliper (0.01 mm). Juveniles and sub-adults were excluded from the analysis according to the complete fusion of cranial sutures (Motokawa et al. 1997, 2003), and by making a histogram of the HBL as an indicator for age identification of small mammals (Li et al. 1989, 1990; Yang 1990).

We calculated the mean and standard deviation of external and skull morphological indices. The pairwise differences between the two species were tested by independent sample t-tests or Mann-Whitney U tests according to results of the Kolmogorov-Smirnov test for their normality of distribution. Principal component analysis (PCA) was used to test the general appropriateness of the groupings supplied by assessment of overall variation in the skull characters. These analyses were performed using SPSS Statistics 24.0 (SPSS, Chicago, IL, USA).

Results

We obtained 1140 bp of mitochondrial DNA sequences from 110 individuals in this study. The ML tree indicated that the specimens we collected were divided into two lineages, one was clustered with the download from GenBank which was distributed in its type locality, i.e., Baoxing of Sichuan Province, and the other was clustered with the download from GenBank which was exclusively distributed in its type locality, i.e., Taiwan Island (Fig. 2). K2P distance of Cytb between these two lineages was 12.3%. Together with the results of blasting on GenBank, a total of 33 specimens of lineage and 77 specimens of lineage collected in this study were genetically identified by Cytb, and their distribution localities plotted in Figure 3. Also, the distribution localities of recently reported in the mainland of China were added to the figure.

Figure 2.

ML tree based on Cytb of genus. Numbers above the branches represent bootstrap support (BS). The blue clade represents and orange clade represents .

Figure 3.

Distributions of and in the mainland of China，Laos and Vietnam. Black and white triangles represent the sampled sites of first presented in this study and in previous studies, respectively. White circles and squares represent the sampled sites of first presented in this study and in previous studies. Black triangles and white circles overlapped indicate sympatry sites.

By investigating our samples of lineage from Baoxing, Sichuan, we found some morphological features correlated with the holotype: the superior articular facets are more angular in dorsal view and the basioccipital region is narrow and ridged particularly anterior to the position of the basioccipital suture in (Fig. 4). On the upper premolar (P4) the protocone is variably positioned relative to the paracone; the posterolingual border of the tooth is not so rounded; and the posterior border of the tooth is deeply concave. The posterobuccal crest of the paracone of the second upper molar (M2) forms a smooth W-shaped loph in unworn dentition (Fig. 5).

Figure 4.

Comparison of crania of (S2576) and (S2566) from Baoxing, Sichuan. Top row from left to right: dorsal views of the skulls of and (S2576 andS2566), ventral views of the skulls in the same order. Lower row: lateral view of skulls and mandibles from top to bottom of and (S2576 and S2566).

Figure 5.

Comparison of teeth of (left: S2576) and (right: S2566) from Baoxing, Sichuan.

A total of 90 adult individuals were screened by age identification including 26 and 64 . The external and two skull measurements (BOC and GPB) were judged as a non-normal distribution by the Kolmogorov-Smirnov test (P<0.05), so we used the Mann-Whitney U Test for interspecific comparisons; for the others with normal distribution (P>0.05) the parametric independent sample t-test was used (Suppl. material 1, Table S3). Descriptive statistics for external and craniodental measurements of the two species and literature measurements (including holotype) are given in Table 1; they were basically consistent with the variation range and limits recorded in the literature except for IOB. was a little larger than in GBL, MPL and BBC. Although there existed significant differences (P<0.05) in some morphological indices between the two species (Table 2), their range of measurements greatly overlapped. In the PCA made on external and skull measurements, three principal components were extracted and captured 70.07% of the total variation. Five indices, GBL, GLS, LUTR, BBC and LM, were the top five with the highest correlations with the first axis (PC1, Table 3). The sample distributions over the scatter plot in coordinate area constructed by first two principal component axes showed a great overlap between the two species in external and skull indices (Fig. 6), indicating that morphological indices cannot accurately identify the two species.

Table 1.

External and cranial measurements of and (in millimeters). Mean ± standard deviation (SD), and range, with number of specimens in parentheses are presented for characters available.

Morphological indices	Crocidura attenuata		Crocidura tanakae
	This study	Jenkins et al. 2009 (including holotype)	This study	Kuroda 1938 (Holotype)	Fang et al. 1997	Motokawa et al. 2001
Total Body Length (TBL)	133.92±6.22	–	132.95±6.84	–	–	–
	120.00–144.00 (26)	–	115.00–148.00 (64)	–	–	–
Head and Body Length (HBL)	77.96±3.54	71.40±6.64	78.90±5.42	66.00	77.90±3.60	73.36±6.67
	71.00–87.00 (26)	63.00–82.00 (10)	70.00–95.00 (64)	–	69.70–86.00	64.5–84.5
Ear Length (EL)	8.32±1.05	–	8.79±1.00	–	8.96±0.39	9.05±0.91
	6.05–10.16 (26)	10.00 (2)	5.97–11.06 (63)	–	8.10–9.80	7.9–10.2
Greatest Length of Skull (GLS)	20.99±0.59	20.55±0.66	20.54±0.61	20.00	20.84±0.42	20.85±0.41
	20.12–22.36 (23)	19.60–21.70 (9)	19.23–21.69 (61)	–	20.03–21.82	19.94–21.31
Cranial base Length (GBL)	20.91±0.57	–	20.39±0.61	–	19.88±0.46	–
	20.04–22.26 (24)	–	18.85–21.68 (61)	–	19.06–20.73	–
Median palatal Length (MPL)	9.07±0.41	–	8.89±0.31	–	–	–
	8.35–10.13 (26)	–	8.14–9.61 (63)	–	–	–
Length of tooth row (LUTR)	9.36±0.26	8.72±0.38	9.16±0.34	9.00	9.29±0.16	9.33±0.19
	8.97–9.96 (26)	8.20–9.50 (11)	8.19–9.91 (63)	–	8.73–9.62	8.87–9.54
Greatest palatal breadth (GPB)	6.25±0.16	6.09±0.22	6.32±0.31	–	6.33±0.16	6.41±0.13
	5.96–6.55 (26)	5.80–6.50 (11)	5.65–6.97 (64)	–	5.98–6.68	6.20–6.58
Breadth of Occipital Condyles (BOC)	5.20±0.19	–	5.13±0.23	–	–	–
	4.78–5.55 (22)	–	4.48–5.64 (60)	–	–	–
Greatest Breath of Braincase (BBC)	9.59±0.27	9.17±0.24	9.26±0.32	–	9.23±0.18	9.22±0.22
	8.98–10.25 (26)	8.70–9.40 (8)	8.57–10.02 (64)	–	8.00–9.66	8.87–9.50
Interorbital Breadth (IOB)	3.97±0.23	4.43±0.18	3.79±0.17	4.50	4.50±0.10	4.62±0.10
	3.61–4.46 (26)	4.10–4.70 (9)	3.43–4.23 (64)	–	4.29–4.68	4.47–4.74
Height of the Braincase (HB)	5.13±0.14	4.91±0.10	5.01±0.15	–	–	–
	4.80–5.42 (26)	4.80–5.10 (8)	4.65–5.35 (64)	–	–	–
Length of mandible (LM)	10.01±0.28	12.76±0.62	9.85±0.33	–	8.38±0.20	–
	9.34–10.67 (26)	11.7–13.9 (11)	9.06–10.52 (64)	–	7.88–8.91	–

Table 2.

Morphological comparisons and significant difference between and in this study. Values in bold show significant differences.

Morphological indices	T test/Mann-Whitney U Test
	F	Sig.	t	Df (N)	Z	P
Total Body Length (TBL)				(90)	-0.687	0.492
Head and Body Length (HBL)				(90)	-0.478	0.633
Ear Length (EL)				(89)	-2.08	0.038
Greatest Length of Skull (GLS)	0.149	0.701	-3.045	82		0.003
Cranial base Length (GBL)	0.133	0.717	-3.541	83		0.001
Median palatal Length (MPL)	1.214	0.274	-2.263	87		0.026
Length of teeth row (LUTR)	1.174	0.282	-2.615	87		0.011
Greatest palatal breadth (GPB)				(90)	-1.077	0.281
Breadth of Occipital Condyles (BOC)				(82)	-1.256	0.209
Greatest Breath of Braincase (BBC)	1.416	0.237	-4.666	88		<0.001
Interorbital Breadth (IOB)	2.492	0.118	-4.143	88		<0.001
Height of the Braincase (HB)	1.131	0.29	-3.61	88		0.001
Length of mandible (LM)	2.766	0.1	2.195	88		0.031

Table 3.

Principal component loadings as performed by analyses of 13 morphological measurements of and .

Variable	Component
	PC1	PC2	PC3
Total Body Length (TBL)	0.691	0.550	0.065
Head and Body Length (HBL)	0.435	0.607	0.458
Ear Length (EL)	0.101	0.706	-0.314
Greatest Length of Skull (GLS)	0.948	-0.088	-0.062
Cranial base Length (GBL)	0.966	-0.096	-0.055
Median palatal Length (MPL)	0.815	0.041	-0.372
Length of teeth row (LUTR)	0.867	-0.139	0.017
Greatest palatal breadth (GPB)	0.566	0.138	0.458
Breadth of Occipital Condyles (BOC)	0.545	-0.358	0.561
Greatest Breath of Braincase (BBC)	0.837	-0.140	-0.022
Interorbital Breadth (IOB)	0.470	0.035	-0.368
Height of the Braincase (HB)	0.590	-0.223	-0.252
Length of mandible (LM)	0.820	-0.126	-0.045
% of total variance explained	49.876	11.087	9.110
Eigenvalue	6.484	1.441	1.184

Figure 6.

Scatter plot of (red) and (blue) sample distribution over PC1 and PC2 axes constructed based on external and skull morphological variables. Different symbols represent different populations. GD: Guangdong, FJ: Fujian, ZJ: Zhejiang, SCBX: Baoxing, Sichuan, HB: Hubei, SC: Sichuan, GX: Guangxi, AH: Anhui, HuN: Hunan, HN: Hainan, CQ: Chongqing, YN: Yunnan, JX: Jiangxi

Among the localities of our field surveys, was recorded at almost all sites investigated (Fig. 3), whereas was only found in the following six provinces: Sichuan Province (Baoxing), Fujian Province (Mount Wuyi), Hubei Province (Shennongjia), Guangdong Province (Nanling), Jiangxi Province (Mount Jinggang), and Zhejiang Province (Jinhua).

Discussion

This study indicates that and are sympatrically distributed not only in continental Indochina (Jenkins et al. 2009, 2013; Bannikova et al. 2011; Abramov et al. 2012) but also in the mainland of China. The distribution of is apparently limited to only two ranges, i.e., Baoxing of Sichuan to Shennongjia of Hubei and Nanling of Guangdong to Jinhua of Zhejiang; the natural range of this species is much smaller than that of which is distributed almost all over the south of mainland China.

Note that the map of (Fig. 1, left) presented by the IUCN is erroneous due to the regular events of species misidentification of in the mainland of China. The IUCN map mistakenly shows the mixed distributions of both and ; the presented distributions of in Taiwan and the Hainan Islands are erroneous for the same reason. For the distribution map of (Fig. 1, right), the range is not definitively established due to the few districts surveyed and information from more recent records has yet not to be included.

Based on morphological features we found among our samples and the results of its comparisons with type materials of and (Jenkins, et al., 2009, 2013), we consider that the specimens of the lineage should be attributed to , and the other lineage to . Wang (2003) divided into three subspecies in China, including the Himalayan subspecies ( Anderson, 1877) distributed in northwestern Yunnan (Gongshan), the South China subspecies ( Milne-Edwards, 1872) distributed in other parts of mainland China and the Taiwan subspecies ( Kuroda, 1938) distributed on Taiwan Island. It is clear that is actually a valid distinct species (Motokawa et al. 2001), but the other two subspecies still need taxonomical validation by detailed analysis to exclude the possibility of misidentification of specimens. Similarly, the same taxonomic challenge exists for Howell, 1926, the subspecies distributed in the Fujian Province (Smith and Xie 2009). All these subspecies are uncertain because the authors may well have wrongly included specimens of mixed with samples.

There are many research reports listing in the mainland of China. For example, Zhang et al. (1987) investigated ( in original paper) as a host animal of epidemic hemorrhagic fever, Wu (2002) reported population density fluctuation in , and many reports on animal diversity and pathogen host studies involved . Gu et al. (2007) reported that epidemiologic surveillance on leptospirosis in the Anhui Province and the first discovery of a pathogenic strain in the renal of , Wu et al. (2008) made a preliminary comparative anatomical study of digestive tracts between and . Because might have been taxonomically misidentified with in these reports, and our present study demonstrates that is much more widely distributed in the mainland of China, the species “” described in these reports may be in fact or at least contains , results of these studies therefore need re-evaluation.