Microhabitat Preferences of the Wild Population of Liua shihi

根据已有的研究，本研究选择巫山巴鲵聚集分布的大巴山区的4个区县作为研究区域开展野外调查工作，包括重庆市的城口县、巫溪县、巫山县以及湖北省的神农架林区(30°45′-32°12′N，108°15′-110°58′E)。

本研究在城口县、巫溪县、巫山县和神农架林区共选择了14条森林溪流作为研究样线，其中7条位于城口县，3条位于巫溪县，2条位于巫山县，2条位于神农架林区(图 1)。样线分布的海拔区间为766~2 043 m，样线平均宽度为4 m，平均长度为300 m。在2023年4月15-19日对城口县、巫溪县和巫山县选定的12条样线进行野外调查，在2023年5月18-19日对湖北省神农架林区的2条样线进行野外调查。根据巫山巴鲵夜间活动更为活跃的特性，调查选择在太阳落山后进行(20：00后)。具体来说，以每一条样线中有巫山巴鲵分布的点为中心，设置1 m × 1 m的分布样方，对照样方的设定是基于分布样方的位置，使用米尺参照随机数生成器产生的角度(上游：0~180°；下游：181°~360°)和距离(0~1 000 cm)来确定位置，对照样点中需没有巫山巴鲵的分布，有巫山巴鲵分布的样点被舍弃。

依据相关文献记载，小鲵科物种的栖息地主要分布在植被茂密、降雨量充沛、气候温和的高海拔地带^[8]。野外调查发现，巫山巴鲵的生活环境主要涵盖沼泽、溪流、水塘以及地下河等生态系统类型(图 2)。因此，对巫山巴鲵的数量和动态进行精确调查存在一定的难度。为了了解小鲵科的种群数量，研究人员常常采用间接估计法，即在个体集中分布的区域通过目测对个体数量进行估算^[12-13]。尽管间接估计法并不能精确获得小鲵科物种的个体数量，但通过这种方法，仍然可以更快速有效地评估这些物种的种群数量状况。本研究中，针对巫山巴鲵数量较少的样点采用精确计数，对于极个别地形复杂、数量较多的样点采用间接估计法，并通过GPS记录每一个样点的地理坐标。

本研究选择并记录了可能对巫山巴鲵分布和个体数量产生影响的14个微生境环境因子，包括13个数值型变量：海拔(Ele)、气温(At)、空气湿度(Rh)、河宽(Rw)、河深(Rd)、水体流速(Fv)、水温(Wt)、水体pH值(WpH)、水体叶绿素a(Chl.a)、水体溶解氧(Do)、水体电导率(Con)、植被覆盖度(Fvc)和水体浊度(Nt)，以及1个分类变量：基质类型(Subt)。具体的测量方法如表 1所示。

通过对14条样线的详细调查，本研究共获得了60个样方的微生境因子数据，其中34个样方有巫山巴鲵分布，分别记录了每个样方内的个体数量，另有26个对照样方没有巫山巴鲵分布。通过比较分布样方和对照样方之间环境因子的差异，并分析环境因子与样方中的个体数量的关系，解析巫山巴鲵的栖息地微生境偏好。对于无序分类变量，利用拟合优度卡方检验分析分布样方和对照样方在基质类型上是否存在显著差异^[14]。然后采用生境喜好系数(Habitat Preference Index，HPI)分析分布样方中有显著差异的基质类型^[15]。根据野外观察到的实际情况，将栖息地基质类型分为了9种，包括沙(＜1 cm)、砾石(1~10 cm)、卵石(10~20 cm)、石块(＞20 cm)、卵石+石块、卵石+沙、砾石+沙、石块+砾石、卵块+砾石。

HPI用以下公式进行计算：

式中：o为观测值；e为期望值。观测值是实际调查中每种基质类型下巫山巴鲵的数量，期望值是假设巫山巴鲵在每种基质类型下均匀分布的数量。HPI的取值范围为-1~1，当-1.0≤h_HPI＜-0.5时，表示回避(E)；-0.5≤h_HPI＜-0.1时，表示不选择(NS)；-0.1≤h_HPI＜0.1时，表示随机选择(RS)；0.1≤h_HPI＜0.5时，表示选择(S)；0.5≤h_HPI≤1.0时，表示偏好(F)。

对于数值型环境变量，先将两组数据中的13个微生境因子进行min-max标准化，以减小极端数值对数据分析的影响，并提高数据精度和可操作性。使用夏皮罗—威尔克检验法(Shapiro-Wilk)检验每种数值变量是否符合正态分布，对分布样方与对照样方的微生境因子进行差异性检验。若微生境因子符合正态分布则使用T检验(T-test)，若不符合正态分布则使用曼—惠特尼秩和检验(Mann-Whitey U)。“****”表示p＜0.000 1；“***”表示p＜0.001；“**”表示p＜0.01；“*”表示0.01＜p＜0.05；“ns”表示p＞0.05。

为降低微生境因子之间的共线性对数据分析造成的影响，还需进行变量间的相关性分析。对符合正态分布的变量采用皮尔逊相关性检验(Pearson)，对不符合正态分布的变量采用斯皮尔曼相关性检验(Spearman)。本研究保留相关性系数|r|＜0.3的变量，并在|r| ≥0.3的变量中保留生态学意义更大的环境变量。最后采用广义线性模型(Generalized Linear Model，GLM)探讨了微生境因子对巫山巴鲵分布的决定作用。GLM无需要求因变量服从正态分布，而因变量可按指数分布任意服从概率分布，因此可直接适用于一般常态线性模型^[16-17]。GLM中数据结构可以存在非线性与非常数方差，所以它的分布模型比较灵活，可以更好地分析生态学中的相互关系。本研究使用赤池信息准则(Akaike Information Criterion，AIC)逐步向后选择模型，即每次进行AIC分析时，拟合度最好的模型具有最小的AIC值，剔除p值最大的环境因子之后，再次对其进行回归分析，将p值最大的变量逐渐剔除，最终获得AIC最小的一个最优模型。

本研究所有的计算和统计分析都在Excel 2016和R 4.2.3软件中进行。

野外调查显示，巫山巴鲵多栖息于灌木茂盛、有水系的周边生境，多静止于岸边或缓慢移动于水中。从拟合优度卡方检验的结果中显示，巫山巴鲵对生境中基质类型的偏好有显著的选择性(表 2)。在对其HPI计算统计后发现，巫山巴鲵倾向于选择栖息在沙基质较多的区域，回避或不选择卵石、石块与砾石较多的区域(表 3)。

Shapiro-Wilk检验结果表明，在13个数值型环境变量中，只有水体pH值不符合正态分布，其余12个都符合正态分布。对水体pH值的两组数据进行T检验，结果表明(表 4)，巫山巴鲵对栖息地中水体pH值的选择没有显著差异；对其余12个环境变量进行Mann-Whitey U检验，结果表明，巫山巴鲵分布样方和对照样方之间的气温、海拔、水体流速、植被覆盖度、空气湿度和河深没有显著差异，而水体叶绿素a、水体溶解氧、水体浊度、河宽、水体电导率和水温存在显著差异。

由于单个变量做相关性分析没有统计学意义，因此将水体pH值纳入其余12个变量，一同做Spearman相关性检验。对全部微生境因子进行筛选，保留|r|＜0.3的变量。对于|r|≥0.3的变量，参考以往的研究^{[4, 8, 13, 15]}，选择对有尾类分布影响更大的变量，最终将海拔、河深、水温、水体pH值、水体叶绿素a、水体电导率和水体浊度7个数值型变量纳入后续广义线性模型的分析(图 3)。

利用广义线性模型将7个微生境因子通过AIC逐步向后选择模型，最小AIC值的模型结果显示，海拔、水温、河深、水体pH值和水体电导率对巫山巴鲵的个体数量有显著影响。参数预测结果表明(表 5)，巫山巴鲵偏好的微生境表现出更低的水温、更浅的河深、较低的海拔、较低的水体pH值和较高的电导率特征，但值得注意的是，电导率(p=0.136)对巫山巴鲵的栖息地选择没有显著影响。在构建最优线性模型时，层次划分的结果表明巫山巴鲵栖息地微生境因子对其分布数量的解释程度各不相同(图 4)。其中，水温是巫山巴鲵栖息数量最重要的预测因子，其解释率达到了48.94%，其次是水体pH值(21.24%)、海拔(19.27%)、河深(7.44%)，解释率最低的因子是水体电导率，只有3.11%。

本研究解析了巫山巴鲵野外种群的栖息地微生境偏好。对于基质类型偏好的研究结果表明，巫山巴鲵多偏向于栖息在环境基质中沙较多的区域，而回避或不选择含有卵石、石块与砾石的区域。以往的研究表明，在森林溪流生态系统中，沙基质较多的区域通常分布有一定数量的水生植物，这些水生植物可提供给巫山巴鲵适合的避难场所。此外，水生植物生长区域往往伴随分布有更多的底栖生物或其他食物来源，因此能够满足巫山巴鲵的生存需求^[18]。野外调查发现，溪流中卵石、石块与砾石较多的区域一般遮蔽物较少，可能会使巫山巴鲵更容易暴露在捕食者(如鸟类)的视野中。同时，这些区域分布的底栖生物，不足以支持巫山巴鲵的捕食需求^[19]。对于两栖动物而言，庇护场所对其生存至关重要，尤其是那些主要依靠皮肤进行气体交换的动物，比如无肺螈科(Plethodontidae)和四趾螈属(Salamandrina)的物种^[20]。同时，由于不同种类两栖动物的体型大小不同，也会有不同基质类型的选择偏好，如镇海棘螈(Echinotriton chinhaiensis)主要选择石块作为繁殖期间的遮蔽物^[15]，而猫儿山小鲵(Hynobius maoershanensis)会倾向于选择腐殖质泥土作为藏身地^[21]。

对水体微生境因子的分析结果表明，巫山巴鲵的分布样方和对照样方相比，水体叶绿素a、水体溶解氧、水体浊度、河深、水体电导率和水温存在显著差异。其中，叶绿素a含量较高的水域，表示水体中可能有更多的浮游植物，即更高的初级生产力，能给巫山巴鲵提供更丰富的食物资源^[22]。由于巫山巴鲵需要依靠皮肤辅助呼气，水中适宜的溶解氧浓度对其生存至关重要，溶解氧浓度较高的水域可以确保巫山巴鲵的呼吸和生存需求。这一研究结果与放归中国大鲵栖息地偏好的结果类似^[4]，表明水体溶解氧可能是驱动水生有尾类分布的重要因素。水体浊度和水体电导率反映了水体中的悬浮物和电解质含量，浊度和电导率较高的水域，可能表明水体中杂质和盐分含量较高。以往的研究表明，这样的水体有利于水生植物的生长^[23]，因此同样有利于巫山巴鲵获取更多的食物。此外，水温对于水生生物的生命活动和代谢有重要影响。广义线性模型结果显示，水温与巫山巴鲵数量表现出显著的负相关关系，因为较高的水温可能影响有尾两栖类的繁殖行为和繁殖成功率，增加卵的死亡率或降低孵化率，从而限制种群的数量增长^[24]。此外，这一结果同样表明巫山巴鲵偏好较低的水温，与东方蝾螈(Cynops orientalis)偏好较高温度的表现相反^[25]。这一结果可能与物种的系统发育关系有关，也在一定程度上反映了物种的热耐受性。而对于气温、海拔、水体流速、植被覆盖度、空气湿度和河深6个因子，巫山巴鲵的分布样方和对照样方没有显著差异，可能是由于研究区域内这些因子存在交互作用，也可能因为巫山巴鲵在研究区域内对这些微生境因子具有较高的耐受性和生态适应性^[4]。但广义线性结果显示，海拔和水体pH值与巫山巴鲵数量也有较为显著的负相关性，较低的海拔可能意味着气候更温暖、湿润，以及更丰富的植物多样性；而较低的pH值可能对巫山巴鲵的某些竞争物种不利。由于本研究考虑的微生境因子有限，未来应纳入更多的微生境因子以更为准确地了解其栖息地偏好。此外，巫山巴鲵在其繁殖期是否有特殊的栖息地需求，也需要在未来的研究中进一步深入探讨。

巫山巴鲵为了更有效地提升自身的存活率，同时降低自身被捕食的风险，更倾向于选择栖息沙基质较多的区域，回避或不选择卵石、石块与砾石较多的区域。本研究还揭示了驱动巫山巴鲵分布的微生境因子，有利于针对性地制定其栖息保护对策。基于本研究的结果，未来在开展巫山巴鲵栖息地保护工作时，不仅需要维持一定的水源，还需要通过生态恢复、植被保护和水体管理等措施控制水温，保持沙基质，并保证一定的初级生产力。同时，要加强水体保护和管理，减少污染物特别是持久性有机污染物(如化肥、农药等)的排放。此外，还需要加强对物种和环境因子的监测调查力度，并加强相关教育宣传和法律保护。