兽研简要总结
欧洲野猫在匈牙利正面临适宜森林栖息地丧失以及与家猫杂交的威胁。这是一个严重问题，因为野猫在防止小型哺乳动物过度繁殖方面发挥着重要的生态作用。本研究旨在评估匈牙利种群中野猫的遗传状况、空间分布以及家猫影响的程度。我们检查了来自该国不同地区的118个个体。结果显示，虽然许多个体属于野猫类群，但相当一部分是自由放养的家猫或表现出混合血统。不同地区由不同的母系谱系主导，反映了复杂的历史，尽管基础遗传工具无法确定这些过去相互作用的确切时间。我们的研究证实，虽然保护监测至关重要，但管理者在做出关于个体分类的关键决策之前需要更先进的遗传工具。这对社会具有高度价值，因为保护野猫对于维持森林生态系统的长期健康平衡至关重要。
1. 引言
在当前前所未有的灭绝危机中，食肉动物，尤其是猫科动物，在哺乳动物中面临最高的灭绝风险[1]。历史上，欧洲野猫曾遭受严重的狩猎和诱捕。这一因素加上栖息地丧失，导致该物种整个分布区的种群数量下降[2,3,4]。主要威胁包括潜在栖息地的丧失[5,6]、人为导致的死亡——尤其是路杀[7,8,9]和偷猎[9]——以及因遗传漂变导致的等位基因变异丧失[3]。一个主要问题是与家猫的杂交[4,7,10,11,12,13,14,15,16,17,18]，据推测这已导致某些种群的局部灭绝。此外，自由放养的家猫可以通过捕食共享的小型哺乳动物猎物[14,15,19]、剥削性和干扰性竞争[20]、干扰以及各种疾病的传播[21]间接影响野猫。作为关键的哺乳动物中型食肉动物，野猫在生态系统功能中发挥着至关重要的作用，特别是在缺乏顶级捕食者的情况下通过自上而下的控制来调节小型哺乳动物种群[22]。然而，自由放养和流浪家猫种群的增加对全球本土野生动物施加了严重的捕食压力[23]。此外，杂交可以从根本上改变这一重要的生态功能，匈牙利野猫、流浪猫及其杂交种之间显著的营养生态位重叠和食性竞争即为证据[14]。
过去几十年的栖息地丧失和破碎化导致种群遗传多样性下降[4]。遗传研究表明，从家猫到野猫的基因流在欧洲各地强度不一，表现出显著的局部差异，这很可能归因于历史或生态因素[17,24]。例如，在意大利、保加利亚、葡萄牙和德国种群中发现了低水平的杂交[11,25,26,27]，而在苏格兰和匈牙利的野猫种群中检测到了高杂交率[4,24,28]。对欧洲野猫种群的研究已识别出五个不同的生物地理类群：北巴尔干、意大利和西西里、德国中部、中欧以及伊比利亚半岛[26,29]。
关于匈牙利野猫分布的信息有限。其种群在地理上可分为三个主要区域：其中最稳定的区域从巴科尼山脉延伸至曾普伦山脉（与北喀尔巴阡山脉相连）；蒂萨河以东地区及蒂萨河洪泛平原森林（与东喀尔巴阡山脉相连）；以及南部与多瑙河洪泛平原相连的迈切克地区。在这些区域，人类土地利用和集约化农业景观构成双重威胁。一方面，密集的道路网络和人类基础设施严重割裂了适宜的森林栖息地，成为阻碍这些核心野猫种群之间自然基因流动的屏障[5,6]。另一方面，相同的景观结构产生了与定居点和农田接壤的广阔林缘网络；这些林缘为家猫提供了侵入先前隔离的野猫栖息地的通道，从而增加了接触和杂交的机会[6]。文献[4]研究了匈牙利种群的杂交模式，得出结论：该地区的样本由可变表型和基因型的组合构成，这可能是长期杂交和基因渗入的结果。根据[30]利用微卫星位点进行的一项研究，在8%的个体中可检测到杂交；然而，使用连锁位点时该值更高，明确证实了匈牙利种群中杂交的持续存在。这种广泛的基因渗入对实际野生动物管理构成了重大挑战。虽然对流浪猫种群的积极调控在匈牙利野猫保护中具有法律和生态上的必要性，但由于深度的基因渗入，传统的形态学鉴定（表型）已不再能可靠地在野外区分纯种野猫与杂交个体[31,32]。因此，保护管理者日益依赖稳健的遗传学工具来精确识别混血个体，并为有针对性的保护和移除工作提供依据[33]。
生物地理学和系统地理学格局是理解分类群进化历史以及确定适当保护管理单元以实现可持续保护的关键因素[16,34]。欧洲野猫的当代遗传结构是更新世气候振荡的直接结果。这些隔离时期促进了遗传漂变和独特线粒体单倍型的固定[35]。迄今为止，仅有少数分子研究调查了欧洲Felis silvestris的系统发育和系统地理学格局。这些研究主要使用信息性线粒体区域（如控制区或NADH脱氢酶亚基）来区分家猫和欧洲野猫谱系，并检测到后者中大致存在一个大陆分支和一个地中海分支[10,27,36,37,38,39]。对线粒体标记（特别是ND5基因）的广泛分析重建了该物种在其整个分布范围内的详尽生物地理历史。已鉴定出五个主要遗传簇，对应于不同的生物地理群组[29]。这些谱系（称为谱系W）的多样化始于约64,200年前，与晚更新世的冷却阶段相吻合[35]。隔离的苏格兰种群、伊比利亚簇和东南欧群组代表了这一多样化中最早的分支[35]。巴尔干半岛是连接小亚细亚和欧洲的重要迁徙走廊，在末次盛冰期（LGM）期间充当避难所[40]。然而，当代巴尔干野猫种群表现出令人惊讶的低线粒体多样性[35]。一个对比鲜明的案例见于撒丁岛野猫（F. s. lybica）。尽管位于欧洲，遗传学和考古学证据将这些猫的起源追溯至2000多年前人类介导的近东野猫扩散[41]。
鉴定“纯种”欧洲野猫与家猫杂交个体依赖于一组经过验证的线粒体多态性。这些标记在法医学和保护遗传学中被优先使用，因为它们能够检测到即使在与野生雄性多代回交之后的母系基因渗入[42]。ND5基因可以说是线粒体基因组中用于区分野生和家养谱系最具诊断价值的部分。系统发育分析已鉴定出两个主要谱系：谱系W（野生）和谱系D（家养/利比亚）。[35]
本研究的目的是描述匈牙利欧洲野猫种群的当前遗传状况，并量化其中与家猫相关血统的程度。为此，我们利用14个常染色体微卫星标记和线粒体ND5基因669 bp片段，评估了野外采集的匈牙利猫的遗传多样性、空间变异以及与家猫相关的血统。我们检验了以下内容：(i) 野外样本是否与家猫参考样本存在分化，并表现出空间结构化的核变异；(ii) 中间STRUCTURE分配结果是否在各亚区间存在差异；(iii) 该标记组合是否能在模拟中区分亲本、F1、F2和回交基因型；(iv) 线粒体单倍型是否表现出地理异质性。我们预期早期世代杂交个体比晚期回交个体更容易被检测到，并将对线粒体谱系的历史解读视为假设，而非直接的时间推断。
2. 材料与方法
2.1. 样本
2016年至2026年间，从推测为野猫的个体采集了毛发样本（n = 118），主要来自道路致死个体。参考家猫毛发样本（n = 30）来自在布达佩斯及其周边兽医诊所就诊的活体家养猫；在常规兽医诊疗或操作过程中，经主人事先知情同意，以非侵入性方式采集了少量带毛囊的毛发（见知情同意声明）。分析前，所有样本均储存于塑料袋中，置于−20 °C。野外采集样本根据匈牙利主要地理区域和该物种全国分布范围的核心区域被划分为地理组。根据其地理位置（补充表S1），野外采集样本被分为以下组：大匈牙利平原（GHP，n = 35）、北外多瑙（TD-N，n = 9）、南外多瑙（TD-S，n = 14），以及北匈牙利山地地区，包括Börzsöny山（NHM-Bö，n = 13）、Cserhát山（NHM-Cs，n = 14）、Mátra山（NHM-M，n = 15）和Bükk山（NHM-Bü，n = 18）（图1）。
2.2. 实验室操作
使用MagCore HF16自动核酸提取仪和MagCore Genomic DNA Tissue Kit（RBC Bioscience Corp.，台北，台湾），按照制造商方案从样本中提取DNA。基因组DNA样本储存于−20 °C。
对于微卫星分析，我们扩增了14个位点（FCA043、FCA023、FCA097、FCA132、FCA096、FCA223、FCA698、FCA126、FCA229、FCA220、FCA090、FCA008、FCA045、F115），这些位点改编自[43]（补充表S2），使用Qiagen Multiplex PCR Kit（QIAGEN GmbH，希尔登，德国）的多重PCR panel进行。反应体系体积为25 μL，包含Multiplex MasterMix（QIAGEN GmbH，希尔登，德国）、60 ng DNA、优化量的10 μM浓度引物以及无核酸酶水。PCR在LifeECO仪器（Bioer Technology Co.，Ltd.，杭州，中国）上进行。所有多重PCR的PCR程序如下：95 °C 15 min；30个循环的94 °C 30 s、58 °C 30 s、72 °C 60 s，最后在60 °C延伸30 min。片段毛细管分析在Applied Biosystems 3500xL Genetic Analyzer仪器上进行，使用LIZ500大小标准（Applied Biosystems，福斯特城，加利福尼亚州，美国）。微卫星片段长度使用Geneious Prime 2024.0.7（Biomatters Ltd.，奥克兰，新西兰）以Two Surrounding Peaks sizing方法进行判定。
从野外样本中，按照[35]的方案扩增并测序了线粒体NADH脱氢酶亚基5基因（ND5；家猫线粒体基因组中核苷酸13,149–13,983；NCBI参考序列NC_001700）。该mtDNA区域特别适合系统发育重建，因为它具有高度鉴别力，包含七个诊断性突变，可用于区分欧洲野猫（F. silvestris）和家猫/非洲猫（F. catus/lybica）谱系，此外，它不含核线粒体片段（numts），且同塑性低[10,36]。ND5分析使用[42]描述的引物进行（正向：5′-TGC CGC CCT ACA AGC AAT-3′，反向：5′-TAA GAG ACG TTT AAT GGA GTT GAT-3′）。反应体系体积为20 μL，包含3 μL稀释至15 ng/μL浓度的DNA、10 μL Phire Animal Tissue PCR Buffer、0.4 μL Phire HS II DNA聚合酶、各0.5 μL稀释至10 μM的引物以及5.6 μL无核酸酶水。PCR在LifeECO仪器（Bioer Technology Co.，Ltd.，杭州，中国）上进行。PCR程序如下：98 °C 5 min，35个循环的98 °C 10 s、60.8 °C 25 s、72 °C 25 s，最后在72 °C延伸5 min。
2.3. 数据分析
2.3.1. 微卫星分析
微卫星数据分析用于分析电泳图并判读等位基因大小。使用MICRO-CHECKER 2.2.3版程序[44]评估无效等位基因、判读错误和大等位基因丢失。为控制与无效等位基因存在相关的伪影，采用在MICRO-CHECKER 2.2.3[44]软件中实现的Brookfield方法[45]。此外，使用FreeNA v1.0程序[46]估算并比较了考虑和不考虑无效等位基因时的FST置信区间。使用GenAlEx v6.51b2[47]中的G统计功能，以9999次置换计算了8个种群的成对FST值。由于FreeNA计算的无效等位基因ENA校正FST值与原始值量级相同，种群对之间的统计显著性（p值）基于在GenAlEx中对原始数据集以9999次置换进行的精确检验确定。我们使用GenAlEx 6.5版[47]估算各 locus 和种群的等位基因频率、等位基因数（Na）、有效等位基因数（Ne）、观测杂合度和期望杂合度值（分别为HO和HE）。
为评估匈牙利不同野猫种群之间以及野猫种群与家猫之间的潜在混群，我们使用STRUCTURE 2.3.4版[48]中实现的贝叶斯聚类方法确定最可能的遗传聚类数，而不预先定义种群。贝叶斯聚类方法和马尔可夫链蒙特卡洛（MCMC）模拟使用混合模型和相关等位基因频率运行。STRUCTURE对每个K值进行10次重复模拟运行，老化期为250,000次迭代，重复750,000次，遗传聚类数（K）从1到7。
使用StructureSelector软件（https://lmme.ac.cn/StructureSelector/，于2026年10月6日访问）[49]确定最佳聚类数（K）。这包括计算数据对数概率的二阶变化率（ΔK），以及应用0.5至0.8之间的阈值评估四种K估计量统计[50]。在使用StructureSelector确定最可能的遗传聚类数（K）后，使用StructureSelector平台内集成的CLUMPAK流程（https://lmme.ac.cn/StructureSelector/，于2026年10月6日访问）[51]对每个K值的10次独立重复运行结果进行对齐和平均。我们应用LargeKGreedy（https://lmme.ac.cn/StructureSelector/，于2026年10月6日访问）算法，以2000个随机输入顺序计算共识个体成员比例q(i)。这些共识q(i)值随后用于生成最终条形图。对于K = 2的解，当野猫聚类成员系数q(i) ≥ 0.75时，野外样本被归类为野猫相关；当q(i) ≤ 0.25时，归类为家猫相关；中间值（0.25 < q(i) < 0.75）则归类为推定混群。该区间以F1基因型预期的q(i) = 0.5值为中心，对应于野猫混群研究中常规应用的范围[27]。
为研究匈牙利野猫种群的空间遗传变异和模式，进行了空间主成分分析（sPCA）[52]。在该分析中，根据STRUCTURE结果被鉴定为家猫的个体（家猫聚类分配q(i) ≥ 0.75）从自由生活样本中排除，仅保留推定混群个体和野猫个体。分析在R软件4.5.0（R Development Core Team，2020）环境中使用adegenet（2.1.11版）和spdep（1.4.2版）包[52]进行，基于14个微卫星 locus 的数据。与传统主成分分析（PCA）不同，sPCA是一种多变量方法，明确将个体的地理坐标与其遗传图谱一起纳入。为计算空间自相关（Moran's I），定义了基于距离的连接网络。为确保每个采样点至少有一个邻居的连续网络，阈值距离根据最近邻居之间的最大距离确定。在构建连接网络时，最小欧几里得距离设为0。最大阈值距离基于个体间最近邻距离的最大值（k = 1邻居标准）动态计算，从而确保空间连续性，并使每个检查个体在分析中至少有一个邻居。sPCA识别出使遗传方差与空间自相关乘积最大化的合成成分，总结了种群的空间遗传结构。为可视化全局模式，第一全局主成分的值以地理尺度沿连续颜色标尺表示。为统计检验距离隔离（IBD），在R 4.5.0统计环境中使用ade4包[54]以999次置换进行了Mantel检验[53]。该检验检查了遗传距离矩阵（欧几里得距离）与地理距离矩阵之间的相关性，以确定关系的强度（r）和统计显著性（p）。
为评估STR标记集在检测家猫与野猫杂交方面的统计效能，我们使用recom-sim.py [55]模拟了亲本基因型及推定混合基因型。作为这些模拟的基线数据，我们选取了完整的实证微卫星图谱（即无缺失位点），这些图谱在初始STRUCTURE分析中被归属于任一亲本组（野猫相关组或家猫相关组），且成员比例q(i) ≥ 0.75。对于K = 2的解决方案，野外样本在野猫聚类的一致成员系数q(i) ≥ 0.75时被归类为野猫相关，在q(i) ≤ 0.25时被归类为家猫相关，中间值（0.25 < q(i) < 0.75）则被归类为推定混合。该区间以F1基因型预期的q(i) = 0.5值为中心，0.75阈值遵循欧洲野猫杂交研究中的既定实践[27]。从自由生活样本中，基于该归属概率有84个个体被归属于野猫聚类；这些个体连同30只家猫被用作参考种群以生成各种基因型类别。我们共模拟了400个基因型，涵盖以下类别：两个亲本类别、F1（亲本×亲本）、F2（F1×F1）、第一代回交（第一代杂种×亲本）和第二代回交（第一代回交后代×亲本）分别与两个亲本世系（每类50个模拟基因型）。随后对这些模拟的亲本及推定混合基因型进行了STRUCTURE分析。与之前运行一致，我们应用了具有相关等位基因频率的混合模型，使用250,000次迭代的预烧期，随后进行750,000次马尔可夫链蒙特卡洛（MCMC）迭代。然而，在此步骤中，我们将测试的遗传聚类数限制为K = 1和K = 2，对每个K值进行10次独立运行。在每个K值内，根据估计的Ln Prob of Data选择10次独立重复中最可能的一次，最高值代表最优运行。
2.3.2. 线粒体分析
118个野外样本的ND5序列在BioEdit 7.7.1.0版本[56]中进行了检查和校正。核苷酸序列通过MEGA12（12.0.14版本）[57]中实现的ClustalW进行比对，并在需要时手动调整。比对结果修剪至669 bp长度（位置13,243-13,911），以便与[35]鉴定的单倍型进行比较，如其补充表S1所述。使用DnaSP 6.12.03版本[58]计算了不同单倍型数目、单倍型多样性和核苷酸多样性以及多态性位点数目。中性检验统计量（Tajima's D和Fu and Li's F*）使用R环境（4.5.0版本；R Development Core Team, 2020）中的pegas（1.4版本）包[59]计算。对于匈牙利序列，最佳核苷酸替换模型通过MEGA12中实现的最大似然分析使用贝叶斯信息准则（BIC）选择，如建议所述[57]。单倍型之间的关系也通过使用PopART 1.7版本[60]中可用的中位数连接法构建单倍型网络进行计算和可视化，设置ε = 10，颠换/转换权重为3:1。
为研究遗传变异的等级分布，使用R环境（4.5.0版本；R Development Core Team, 2020）中的poppr（2.9.8版本）包[61]计算了分子方差分析（AMOVA）。分析运行10,000次置换以评估每个参数的显著性。对于等级分析，种群按地区分组：TD-S和TD-N归入Transdanubia，NHM-Bö、NHM-Cs、NHM-M和NHM-Bü归入North Hungarian Mountains，GHP归入Great Hungarian Plain。
3. 结果
本研究中调查的野外野猫个体的详细地理归属、核微卫星分析得出的遗传混合比例（q(i)值）以及鉴定的线粒体单倍型和单倍群总结于补充表S1中。
3.1. 微卫星
对于家猫和野猫，所有研究的位点均具有多态性。在野猫组中，单位点等位基因数最低为3个（NHM-M组的FCA220位点），最高为14个（GHP组的FCA223和F115位点）；对全部14个位点取平均，每组平均等位基因数为7.07至10.79。极端值出现在不同位点——仅在NHM-M组内，FCA220携带3个等位基因而FCA096携带13个——因此单位点值和平均值并不矛盾。平均观测杂合度（HO）为0.73–0.82（在北匈牙利山脉—伯尔热尼组和外多瑙南部组），而平均期望杂合度（HE）为0.76–0.82（在北匈牙利山脉—比克组和大匈牙利平原组）。观测杂合度（HO）仅在南部外多瑙（TD-S）和比克组（NHM-Bü）高于期望杂合度（HE），而在所有其余野猫组中HE均高于HO。有效等位基因数为4.56至5.88（补充表S3）。单位点上HO与HE之间的偏离呈双向性——54%的位点×组组合显示HO < HE，46%显示HO > HE——对14个位点取平均，两个统计量在每组中均非常接近，野猫组平均杂合子缺陷为3.5%。可重复的缺陷仅限于少数位点，尤其是FCA220，这与将家猫开发的标记应用于野猫时所预期的无效等位基因一致；这些位点用MICRO-CHECKER进行了筛查，并用FreeNA中实施的ENA校正重复了群体比较，所得FST值与未校正估计值处于同一数量级。将空间结构化的组视为单一群体（瓦伦德效应，与区域内群体间的显著分化一致，ΦSC = 0.083，p = 0.0226）也会导致所观察到的小正值。
STRUCTURE程序显示两个遗传簇时平均对数似然值最高，对数似然值的二阶变化率也以两个遗传簇时最高，K = 2（图2A）。在这种情况下，家猫参考样本聚为一组，而各个野猫组聚为另一组。在地理定义的野猫组内，18个个体在遗传上被鉴定为家猫簇，16个个体被鉴定为推定混合个体。Puechmaille方法的四项监督统计量表明在所有阈值下均存在三个簇，K = 3（图2B）。在这种情况下，家猫和野猫形成不同的遗传簇；然而，在野猫组内观察到的内部遗传结构（K = 3）与地理位置不相关。由于这些内部野猫遗传谱系在研究区域内未形成地理隔离的群体，所有后续群体遗传分析均在K = 2模型下进行，假设野猫代表单一多态性群体。在地理定义的野猫组内，该模型下也鉴定出家猫（n = 18）和推定混合个体（n = 10）（图2B）。
基于Structure聚类分析的结果，NHM-Bü亚区在所调查的地理组中推定混合个体和家猫个体数量最低（两类均为n = 0）。相反，推定混合个体数量最高的是GHP群体（n = 7；20%），该群体家猫出现率也最高（n = 6；17.14%）。检查其他区域，来自TD-N区域的9个样本中，2个（22.22%）在遗传上被归为野猫，4个（44.44%）为推定混合个体，3个（33.33%）为家猫。来自TD-S区域的14只动物绝大多数为野猫（n = 12；85.71%），同时仅鉴定出1个推定混合个体（7.14%）和1只家猫（7.14%）。在来自NHM-Bö组的13个分析样本中，9个（69.23%）被归类为野猫，而推定混合个体和家猫各2个（15.38%）。在来自NHM-Cs组的14个样本中，11个（78.57%）为野猫，1个为推定混合个体（7.14%），2个为家猫（14.29%）。最后，在NHM-M区域的15个个体中，10个（66.67%）代表野猫谱系，另外检测到1个推定混合个体（6.67%）和4只家猫（26.67%）（图3）。
基于空间主成分分析（sPCA）的结果，匈牙利野猫种群的遗传变异并非在空间上随机分布。在图4中可以直观地识别出空间遗传结构。空间上接近的个体（由颜色深浅相近的点表示）表现出较高的遗传相似性，而位于该国较远地区的种群之间则存在明显的遗传分化。这种直观识别的空间结构得到了Mantel检验结果的统计学证实。遗传距离与地理距离之间检测到显著正相关（r = 0.118，p = 0.005）。这种关联在统计学上显著但较弱，与采样种群内距离隔离（IBD）的模式一致，但并不能确证该模式的存在。
所研究种群间的遗传距离通过经无效等位基因校正（FreeNA）的成对FST值进行估计（表1）。结果显示，家猫种群与所有其他研究群体存在强烈且高度显著的遗传分化（FST = 0.034–0.139；p < 0.01和p < 0.001）。整个数据集中分化程度最高的是家猫与NHM-Bü种群之间（FST = 0.139，p < 0.001），很可能是因为在该群体中未发现家猫或推定混合个体。非家猫群体之间的遗传结构较为中等但界定清晰。NHM-Bü种群与其他群体表现出统计学显著的中等分化（FST = 0.026–0.062；p < 0.05–0.001）。相比之下，在某些地理或分类种群对之间——特别是NHM-Bö、NHM-Cs和NHM-M群体之间——成对FST值与零无显著差异（FST ≤ 0.017；p > 0.05；负值估计按零处理）。因此，这些群体之间未检测到显著分化，这与它们之间的基因流一致；然而，鉴于群体规模小且不均衡，该结果不应被解读为基因流无屏障的正面证据。
非对角线元素表示指定种群群体之间的遗传距离。星号表示统计显著性水平：*** p < 0.001，** p < 0.01，* p < 0.05，而“ns”表示无显著分化（p > 0.05）。GHP（大匈牙利平原）、TD-N（北外多瑙）、TD-S（南外多瑙）、NHM-Bö（北匈牙利山地—Börzsöny）、NHM-Cs（北匈牙利山地—Cserhát）、NHM-M（北匈牙利山地—Mátra）和NHM-Bü（北匈牙利山地—Bükk）。
使用先前应用的分配概率q(i) ≥ 0.75，生成的野猫基因型100%分配到一个聚类，家猫基因型100%分配到另一个聚类。在生成的F1和F2中，分别有100%和90%的基因型被鉴定为推定混合个体；F1基因型均未分配到亲本聚类，而其余10%的F2基因型被分配到家猫聚类。在生成的第二次回交中，76%和84%分配到相应的亲本聚类，而24%和16%未能分配。第三次回交100%分配到适当的亲本群体，没有生成的基因型未被分配。关于生成基因型的Structure q(i)值的详细信息见补充表S4和S5。
3.2. 系统发育分析
通过分析669 bp的序列片段，我们在118个样本中鉴定出八种单倍型。比对序列包含16个多态位点，对应16个突变；在不含重组的无限位点模型下，该数量的分离位点最多可支持17种不同的单倍型，此处观察到其中八种。分析数据集的总单倍型多样性（Hd）为0.683 ± 0.036，核苷酸多样性（π）为0.00704 ± 0.00062。个体间平均核苷酸差异数（k）为4.712。整体数据集的人口统计学统计量未显示显著偏离突变-漂变平衡，Tajima's D值为1.5719（p = 0.1160），Fu和Li's F*值为0.6815（p = 0.4617）（表2）。
鉴定出的八种单倍型属于四个主要单倍群：亚分支Ia（n = 60；50.85%）以优势单倍型dw4为代表，其次是单倍群III（n = 25；21.19%，主要由单倍型dw1驱动，n = 22；18.64%）、单倍群V（n = 15；12.71%；w2）、单倍群IV（n = 12；10.17%；w1）和亚分支Ib（n = 6；5.08%；dw3），而单倍型dw5仅在单个样本中发现（0.85%）。鉴定出两种新单倍型（Hap1和Hap7），均属于单倍群III，各由一个样本代表（0.85%）。
当数据集被划分为七个地理亚区时，鉴定出的单倍型数量范围为三个（TD-S、NHM-Cs和NHM-Bü）至五个（GHP）（图5）。单倍型dw4在NHM种群（NHM-Bö、NHM-Cs、NHM-M和NHM-Bü）中占优势，而dw1在GHP组中最常见，dw1和dw4在TD-N中同样普遍。相比之下，TD-S种群以单倍型w1占优势为特征。遗传变异性在TD-N和GHP亚区最高，无论是单倍型多样性还是核苷酸多样性，而最低的多样性值在NHM-Cs和NHM-Bö组中观察到（表2）。人口统计学统计量（Tajima's D和Fu and Li's F*）未显示任何亚区显著偏离平衡（p > 0.05），尽管NHM-Bö组显示出接近显著性的负趋势（Tajima's D = −1.7203，p = 0.0854；Fu and Li's F* = −1.9014，p = 0.0699）。
3.3. 匈牙利野猫及推定混合个体的系统地理学分析
重建的MJ网络包含从分析数据集中鉴定出的八种单倍型（dw1；dw3；dw4；dw5；w1；w2；Hap1；Hap7），揭示了明显的遗传结构。网络拓扑结构围绕两个高频且占优势的单倍型构建：dw4和dw1。这两个中心变异体通过一个中间的低频单倍型（dw3）相隔总共三个突变步骤，dw3与dw4相隔一个突变，与dw1相隔两个突变。节点dw1形成星状辐射的中心：低频的、可能为私有或稀有的单倍型（Hap1、Hap7、dw5）通过一个突变步骤从该中心变异体分化出来。相反，w1和w2代表一个遗传上更远的分化分支。w1与中心dw1变异体相隔七个突变步骤，随后还有三个额外突变导致单倍型w2。颜色编码的地理分布突出显示，最频繁的单倍型（特别是dw4和dw1）并非严格局限于单一地理区域，而是在多个亚区之间共享（例如，大匈牙利平原/GHP和各个北匈牙利山脉/NHM组）。这种分布模式支持种群之间存在遗传连通性和历史基因流的观点（图6）。
为研究遗传变异的层级分布，进行了分子方差分析（AMOVA），揭示绝大部分变异分布在样本内（σ2 = 3.899；占总变异的77.88%；ΦST = 0.221；p = 0.0001）。区域内种群间的分化占方差的一小部分但具有统计学显著性（σ2 = 0.354；7.08%；ΦSC = 0.083；p = 0.0226）。相比之下，三个地理大区之间的分化解释了15.04%的变异，且不具有统计学显著性（σ2 = 0.753；ΦCT = 0.150；p = 0.1545）（表3）。
4. 讨论
用于匈牙利野猫种群遗传分析的微卫星位点——其中大多数最初是为家猫开发的[43]——在估计种群遗传指标方面表现良好。匈牙利样本中测得的平均等位基因数（NA：7.07–10.79）和期望杂合度（HE：0.76–0.82）较高。这种等位基因多样性超过了其他几个欧洲种群报道的值，如德国（HE：0.65）和意大利（HE：0.69）[62,63]。匈牙利的结果与法国（NA：7.0；HE：0.73）和克罗地亚（NA：7.2；HE：0.72）观察到的值最为相当[64,65]。然而，必须注意的是，从自然栖息地采集的野猫样本中观察到的看似较高的多样性指数往往并不反映种群“健康”，而是提示存在来自家猫的历史或近期基因渗入（基因流），如[33]在瑞士种群中所证明的那样。
尽管STRUCTURE分析（在K = 3时）在自由生活野猫样本中识别出了聚类，但这些聚类未表现出任何明显的空间模式。这些结果表明某些区域（如Börzsöny、Cserhát、Mátra）之间的遗传分化有限，提示区域分化较弱、基因流 substantial 且存在持续混合。东欧——特别是匈牙利——野猫种群在遗传上高度混合和整合；因此，无法像苏格兰等种群那样被记录的那样分化为不同的生物地理单元。这可能是长期杂交和基因渗入的结果[4]。虽然在其他欧洲种群中观察到了由山脉或主要公路引起的明显遗传不连续性[62,66]，但在匈牙利的研究区域之间不存在这样的显著自然或人为基因流屏障。
相比之下，sPCA揭示了高分辨率结构，因为该方法纳入了样本的地理坐标。在移除已识别的家猫个体后，Mantel检验表明遗传距离与地理距离之间存在显著但微弱的正相关（r = 0.118，p = 0.005）。这表明存在距离隔离。匈牙利种群的特征不是强烈的遗传断裂，而是逐渐的、连续的过渡（渐变群），这可能提示栖息地连通性良好和该物种的高移动性[66]。
基于这些结果可以认为，匈牙利自然种群中存在野猫与家猫的杂交。在我们当前的研究中，基于STRUCTURE分析检测到13.56%的中间分配率，而野外鉴定的自由放养家猫又占样本的15.25%（使非野猫合计比例达到28.81%）。作为方法学局限，应当指出，用于模拟的亲本野猫参考个体同样是基于STRUCTURE q(i) > 0.75阈值选出的，而非经全基因组基因组学方法独立验证的纯野猫参考。这与欧洲其他多个种群的研究结果一致；然而，杂交程度似乎在不同地区差异显著。德国的杂交率相对较低，估计约为3%[27,32]，而其他国家记录了提示中度杂交的较高水平，如意大利（3–11%）、瑞士（12.5–15%）和伊比利亚半岛（14–21%）[11,18,32,37,62,67]。参考文献[64,65]报告法国的杂交率约为29%，巴尔干各种群记录的数字为13%至52%。最极端的杂交率见于苏格兰，数值可达80–90%[12,24]。在匈牙利，2000年代中期进行的初步研究使用遗传学方法估计杂交率在25%至31%之间[4,30]。
在18世纪和19世纪，欧洲野猫遭到密集狩猎，导致种群急剧下降、局部灭绝和破碎化[68]。当种群密度达到临界水平时，同种个体的稀少增加了与无处不在的家猫交配的可能性[4]。森林砍伐和栖息地破碎化迫使野猫进入小而孤立的栖息地[68]。这导致遗传多样性下降，并增加了杂交风险，尤其是在栖息地边界处[30]。农业区和森林斑块构成的镶嵌景观（农林复合生态系统）促进了家猫与野猫的相互作用。连续森林为野猫提供庇护，而居民点和农业用地则成为家猫和杂交个体的中心；在这些接触区内，杂交变得更加频繁[17,24]。在匈牙利野猫中观察到的广泛杂交，历史上早已由迫害导致的栖息地丧失和种群下降所注定。这场人口危机，加之匈牙利特有的镶嵌景观结构，使得野猫与家猫在森林边缘的接触增加。结合我们当前的遗传学发现，这些因素证实了基因渗入的古老性和持续性，而这基本上威胁着当地种群的遗传完整性[4,30]。
除遗传纯度丧失外，这种持续的杂交还可能产生严重的生态后果。作为中型食肉动物，野猫在森林生态系统中对小型哺乳动物种群提供重要的自上而下控制[22]。广泛的杂交可能改变表型特征、狩猎行为和栖息地利用，可能使杂交个体的生态位向野化家猫的生态位靠拢，从而破坏其原生生态系统功能[14]。因此，迫切需要严格且有针对性保护管理。由于我们的空间分析（sPCA）揭示的是逐渐过渡的渐变群而非急剧的遗传断裂，管理策略必须聚焦于农林接触区。为减轻进一步的基因渗入，欧洲各地区的保护策略通常侧重于在核心森林栖息地周围建立缓冲区[6,32]。在这些接触区内，野生动物管理实践常涉及对自由放养和野化家猫的调控或清除[32]。然而，由于形态学鉴定在高度杂交种群中已被证明不可靠[31]，若干国家的此类清除工作正日益得到快速遗传筛查的支持，以避免误杀纯野猫[33]。此外，国际保护框架强调，针对农村猫主人的公众意识宣传活动，宣传负责任养宠（如绝育和室内饲养），是旨在减少家猫基因不断流入野猫基因库的综合管理计划的重要组成部分[12,27]。
对匈牙利野猫种群线粒体DNA的分析揭示了八种不同单倍型的存在，这些单倍型构成了以两个优势变异体dw4和dw1为中心的结构化中接网络。在更广泛的欧洲系统地理学背景下评估，这些单倍型明显分为两个深度分化的进化谱系：谱系D，包含家猫和非洲野猫的起源；谱系W，代表归属于野猫的祖先成分[10,35,36]。共享（dw）单倍型在匈牙利各子区域的广泛分布最初似乎证实了家猫在全国范围内的广泛基因渗入。然而，详细的系统发育评估揭示了更为复杂的进化历史。匈牙利种群中最显著的系统地理格局之一是dw4单倍型在北匈牙利山地地区（NHM-Bö、NHM-Cs、NHM-M和NHM-Bü）的压倒性优势。该变异体属于谱系D的亚分支Ia，根据进化模型，其分化时间约为37,700年前[35]。由于这一分化时间远早于目前已知的猫驯化考古学和遗传学证据[10,69,70]，dw4在北部栖息地的高频率并不一定表明与现代家猫的近期基因渗入。相反，一种可能的解释是晚更新世历史谱系保留和古代基因渗入的遗产，当时向黎凡特和安纳托利亚迁移的欧洲野猫与当地Felis lybica种群杂交[36,71,72]。这种古代基因渗入可能反映了在北匈牙利山地地区观察到的线粒体-核基因不一致现象，该地区的个体基于核微卫星看似遗传纯正，却携带共享的线粒体模式[26]。相比之下，主要发现于大匈牙利平原（GHP）的dw1单倍型揭示了完全不同的进化故事。基于人口统计学模拟，dw1和dw3等单倍型在约1000年前近期分化，可能提示近期由人类介导的源自Felis catus的基因渗入[35]。这种近期基因渗入与大匈牙利平原核STR分析检测到的高比例混血个体完全一致，反映了在高度破碎化的农林镶嵌景观中强烈的人为压力和与家猫的频繁互动[4,6]。欧洲野猫（谱系W）在匈牙利种群中的祖先进化遗产主要由w1和w2单倍型保存。值得注意的是，w1单倍型在南外多瑙河地区（TD-S）尤为常见，该地区基于核结果也表现出高比例的归属野猫个体。根据大陆系统地理学框架，w1变异体属于分支IV，对应于东南欧空间集群[35]。这一归属于野猫的谱系在匈牙利南部的优势可能反映了更新世气候振荡后源自巴尔干或地中海以外冰川避难所的冰后期重新定殖浪潮[4,29,73]。这些遗传格局的复杂而连续的性质得到了分子方差分析（AMOVA）的进一步支持。遗传变异的层次划分显示，绝大多数变异（77.88%）存在于个体样本内，而主要地理区域间的遗传分化（15.04%）不具有统计学显著性。缺乏明显的大尺度区域隔离与中接网络的拓扑结构相关，后者清楚地表明优势单倍型在地理上遥远的子区域之间大量共享。总体而言，这些结果表明匈牙利野猫种群并未以长期地理隔离为特征；其遗传结构一直受到该物种高流动性、历史重新定殖浪潮、古代分类单元间基因流以及近期人为杂交在不同地区强度差异的持续塑造[35,66]。
5. 结论
本研究对匈牙利欧洲野猫种群的遗传多样性、空间结构以及杂交带来的生态挑战进行了全面概述。在北匈牙利山脉观察到的与历史谱系保留和古老基因渗入相符的模式（推测源于过去与Felis lybica谱系的混合），以及在大匈牙利平原检测到的与家猫的近期杂交，充分反映了该种群复杂的进化和种群历史。空间分析（sPCA）未揭示匈牙利采样栖息地之间存在明显的遗传屏障；该种群以连续的空间过渡和距离隔离为特征，与景观尺度上的空间基因流相一致。然而，广泛的杂交严重威胁野猫的遗传完整性及其在森林生态系统中调控小型哺乳动物种群的生态作用。由于高度混合，基于形态特征的鉴定已不再可靠。因此，未来的保护策略必须侧重于在关键农林接触区建立缓冲区，通过快速可靠的遗传筛查方法支持流浪家猫的管控，重新审视当前对自然出现的推定混合个体的法律规定，并开展关于负责任宠物饲养（如绝育、室内饲养）的公众意识宣传活动，这些对于阻止家猫基因的持续流入至关重要。微卫星位点和线粒体ND5基因的联合分析表明，采样种群在所检测的标记上保留了较高的遗传多样性，且与历史和近期基因渗入的同时存在相符。鉴于样本量和标记组合的分辨率，这些模式应被视为指示性的而非结论性的。特别是，样本量较小的子区域（如TD-N，n=9）的多样性估计存在较大的抽样误差，应谨慎解读，而不应作为直接区域比较的依据。
查看英文原文▾
Simple Summary
The European wildcat is threatened by the loss of suitable forest habitats and hybridisation with domestic cats in Hungary. This is a severe issue because wildcats play an important ecological role in preventing the overpopulation of small mammals. Our research aimed to assess the genetic status, spatial distribution, and extent of domestic cat influence within the Hungarian population. We examined 118 individuals from various regions of the country. Our results revealed that while many individuals belong to the wildcat group, notable proportion are free-ranging domestic cats or show mixed ancestry. Different maternal lineages dominate in different regions, reflecting a complex history, although basic genetic tools cannot pinpoint the exact timing of these past interactions. Our study confirms that while conservation monitoring is essential, managers need more advanced genetic tools before making critical decisions about individual classification. This is highly valuable to society, as protecting wildcats is fundamental to maintaining the long-term, healthy balance of forest ecosystems.
1. Introduction
In the current unprecedented extinction crisis, carnivores, and felids in particular, are exposed to the highest extinction risk among mammals [1]. Historically, the European wildcat has been subjected to severe hunting and trapping. This factor, combined with habitat loss, has resulted in population declines throughout the species’ entire distribution [2,3,4]. The primary threats include the loss of potential habitats [5,6], human-induced mortality—particularly roadkills [7,8,9] and poaching [9]—as well as the loss of allelic variants due to genetic drift [3]. A major issue is hybridisation with the domestic cat [4,7,10,11,12,13,14,15,16,17,18], which is hypothesized to have already led to the local extinction of certain populations. Furthermore, free-ranging domestic cats can affect wildcats indirectly through predation on shared small-mammal prey [14,15,19], exploitative and interference competition [20], disturbance, and the transmission of various diseases [21]. As a key mammalian mesocarnivore, the wildcat plays a crucial role in ecosystem functioning, particularly in regulating small mammal populations through top-down control in the absence of apex predators [22]. However, the increasing population of free-ranging and feral domestic cats exerts severe predation pressure on native wildlife globally [23]. Furthermore, hybridisation can fundamentally alter this essential ecological function, as evidenced by significant trophic niche overlaps and dietary competition between wildcats, feral cats, and their hybrids in Hungary [14].
Habitat loss and fragmentation over the past decades have led to a decrease in the genetic diversity of populations [4]. Genetic studies have revealed that gene flow from domestic cats to wildcats varies in intensity across Europe, exhibiting significant local differences that can likely be attributed to historical or ecological factors [17,24]. For instance, low levels of hybridisation have been found in Italian, Bulgarian, Portuguese, and German populations [11,25,26,27], whereas high hybridisation rates have been detected in the wildcat populations of Scotland and Hungary [4,24,28]. Examining the European wildcat population, five distinct biogeographical groups have been identified: the Northern Balkans, Italy and Sicily, Central Germany, Central Europe, and the Iberian Peninsula [26,29].
Information regarding the distribution of the wildcat in Hungary is limited. The population can be geographically divided into three main areas: the most stable of these spans from the Bakony Mountains to the Zemplén Mountains (connecting to the Northern Carpathians); the Transtisza region and the floodplain forests of the Tisza River (connecting to the Eastern Carpathians); and in the south, the Mecsek region connected to the Danube floodplain. In these regions, human land use and intensive agricultural landscapes pose a dual threat. On the one hand, a dense road network and human infrastructure heavily fragment suitable forest habitats, acting as barriers that hinder natural gene flow between these core wildcat populations [5,6]. On the other hand, the same landscape structure produces an extensive network of forest edges bordering settlements and farmland; these edges provide routes along which domestic cats can penetrate into previously isolated wildcat habitats, thereby increasing the opportunity for contact and hybridisation [6]. Investigating the hybridisation patterns of the Hungarian population, ref. [4] concluded that the samples from this region comprise a combination of variable phenotypes and genotypes, which may be the result of long-term hybridisation and introgression. Based on a study utilizing microsatellite loci by [30], hybridisation was detectable in 8% of the individuals; however, using linked loci, this value was higher, clearly confirming the ongoing presence of hybridisation in the Hungarian population. This extensive introgression presents a significant challenge for practical wildlife management. While the active regulation of feral cat populations is legally and ecologically essential for wildcat conservation in Hungary, traditional morphological identification (phenotype) is no longer reliable for distinguishing pure wildcats from hybrids in the field due to advanced introgression [31,32]. Consequently, conservation managers are increasingly dependent on robust genetic tools to precisely identify admixed individuals and to inform targeted conservation and removal efforts [33].
Biogeographic and phylogeographic patterns are key factors in understanding the evolutionary histories of taxa and in identifying appropriate conservation management units for sustainable preservation [16,34]. The contemporary genetic structure of the European wildcat is a direct consequence of the climatic oscillations of the Pleistocene epoch. These periods of isolation promoted genetic drift and the fixation of unique mitochondrial haplotypes [35]. To date, only a few molecular studies have investigated the phylogenetic and phylogeographic patterns of Felis silvestris in Europe. These investigations mainly used informative mitochondrial regions (such as the Control Region or the NADH dehydrogenase subunits) to distinguish domestic cat and European wildcat lineages, and detected the rough presence of a continental and a Mediterranean clade within the latter [10,27,36,37,38,39]. Extensive analysis of mitochondrial markers, particularly the ND5 gene, has reconstructed an exhaustive biogeographic history of the species across its entire distribution range. Five main genetic clusters have been identified, corresponding to distinct biogeographic groups [29]. The diversification of these lineages (referred to as Lineage W) began approximately 64,200 years ago, coinciding with the cooling phases of the late Pleistocene [35]. The isolated Scottish population, the Iberian cluster, and the South-Eastern European group represent the earliest branches of this diversification [35]. The Balkan Peninsula has served as a critical migration corridor between Asia Minor and Europe, acting as a refugium during the Last Glacial Maximum (LGM) [40]. However, contemporary Balkan wildcat populations exhibit surprisingly low mitochondrial diversity [35]. A contrasting case is found in the Sardinian wildcat (F. s. lybica). While located in Europe, genetic and archeological evidence traces the origins of these cats back to human-mediated dispersal of Near Eastern wildcats over 2000 years ago [41].
The identification of “pure” European wildcats versus domestic hybrids relies on a set of vetted mitochondrial polymorphisms. These markers are prioritized in forensic and conservation genetics because of their ability to detect maternal introgression even after multiple generations of backcrossing with wild males [42]. The ND5 gene is arguably the most diagnostic portion of the mitogenome for distinguishing between wild and domestic lineages. Phylogenetic analyses have identified two primary lineages: Lineage W (Wild) and Lineage D (Domestic/Lybica) [35].
The aim of this study was to characterize the current genetic status of the Hungarian European wildcat population and to quantify the extent of domestic cat associated ancestry within it. To this end, we assessed genetic diversity, spatial variation, and domestic cat associated ancestry in field collected Hungarian cats using 14 autosomal microsatellites and a 669 bp fragment of mitochondrial ND5. We tested whether (i) field samples were differentiated from domestic cat references and showed spatially structured nuclear variation; (ii) intermediate STRUCTURE assignments varied among sub-regions; (iii) the marker panel could distinguish parental, F1, F2, and backcross genotypes in simulations; and (iv) mitochondrial haplotypes exhibited geographic heterogeneity. We expected early-generation hybrids to be more readily detected than advanced backcrosses and treated historical interpretations of mitochondrial lineages as hypotheses rather than direct temporal inferences.
2. Materials and Methods
2.1. Samples
Hair samples (n = 118) were collected between 2016 and 2026 from putative wildcats, predominantly from road-killed individuals. Reference domestic cat hair samples (n = 30) were obtained from live, owned domestic cats treated at veterinary clinics in and around Budapest; a small number of hairs with follicles was taken non-invasively during routine veterinary care or handling, with the prior informed consent of the owners (see Informed Consent Statement). Prior to analysis, all samples were stored in plastic bags at −20 °C. The field-collected samples were assigned to geographic groups based on the major geographical regions of Hungary and the core areas of the species’ national distribution range. Based on their geographic locations (Supplementary Table S1), the field-collected samples were categorized into the following groups: Great Hungarian Plain (GHP, n = 35), Northern Transdanubia (TD-N, n = 9), Southern Transdanubia (TD-S, n = 14), as well as the areas of the North Hungarian Mountains, Börzsöny (NHM-Bö, n = 13), Cserhát (NHM-Cs, n = 14), Mátra (NHM-M, n = 15), and Bükk (NHM-Bü, n = 18) Mountains (Figure 1).
2.2. Laboratory Procedures
DNA from samples was extracted using the MagCore HF16 Automated Nucleic Extractor with the MagCore Genomic DNA Tissue Kit (RBC Bioscience Corp., Taipei, Taiwan), according to the manufacturer’s protocol. The genomic DNA samples were stored at −20 °C.
For microsatellite analysis we amplified 14 loci (FCA043, FCA023, FCA097, FCA132, FCA096, FCA223, FCA698, FCA126, FCA229, FCA220, FCA090, FCA008, FCA045, F115) adapted from [43] (Supplementary Table S2), using multiplex PCR panels with the Qiagen Multiplex PCR Kit (QIAGEN GmbH, Hilden, Germany). The reaction mixture volume was 25 μL, containing Multiplex MasterMix (QIAGEN GmbH, Hilden, Germany), 60 ng of DNA, optimized amounts of primers at 10 μM concentration, and nuclease-free water. PCR was conducted in a LifeECO instrument (Bioer Technology Co., Ltd., Hangzhou, China). The PCR profile for all multiplexes was the following: 95 °C for 15 min; 30 cycles of 94 °C for 30 s, 58 °C for 30 s, 72 °C for 60 s, and final extension at 60 °C for 30 min. Fragment capillary analysis was conducted on Applied Biosystems 3500xL Genetic Analyzer instruments with LIZ500 size standard (Applied Biosystems, Foster City, CA, USA). The lengths of microsatellite fragments were described using Geneious Prime 2024.0.7 (Biomatters Ltd., Auckland, New Zealand) with Two Surrounding Peaks sizing method.
From the field samples, the mitochondrial NADH dehydrogenase subunit 5 gene (ND5; nucleotides 13,149–13,983 in the domestic cat mitochondrial genome; NCBI reference sequence NC_001700) was amplified and sequenced following the protocol of [35]. This mtDNA region is particularly suitable for phylogenetic reconstructions because it is highly discriminatory, containing seven diagnostic mutations that are useful for distinguishing between the European wildcat (F. silvestris) and domestic cat/African cat (F. catus/lybica) lineages, and furthermore, it does not contain nuclear mitochondrial segments (numts) and has low homoplasy [10,36]. Analysis of ND5 was done using primers described by [42] (Forward: 5′-TGC CGC CCT ACA AGC AAT-3′, Reverse: 5′-TAA GAG ACG TTT AAT GGA GTT GAT-3′). The reaction mixture volume was 20 μL, it contained 3 μL of DNA diluted to the concentration of 15 ng/μL, 10 μL of Phire Animal Tissue PCR Buffer, 0.4 μL of Phire HS II DNA polymerase enzyme, 0.5 μL of each primer diluted to 10 μM and 5.6 μL nuclease free water. PCR was conducted in a LifeECO instrument (Bioer Technology Co., Ltd., Hangzhou, China). The PCR profile was the following: 98 °C for 5 min, 35 cycles of 98 °C for 10 s, 60.8 °C for 25 s, 72 °C for 25 s, and final extension at 72 °C for 5 min.
2.3. Data Analysis
2.3.1. Microsatellite Analysis
Data analysis for microsatellites was used to analyze electropherograms and score allele sizes. The program MICRO-CHECKER ver. 2.2.3 [44] was used to assess null alleles, scoring errors and large allele dropout. To control for artefacts associated with the presence of null alleles, the Brookfield method [45] as implemented in the MICRO-CHECKER 2.2.3 [44] software was used. Moreover, estimates of FST confidence intervals considering and without considering null alleles were estimated and compared using the program FreeNA v1.0 [46]. Pairwise FST values for the eight populations were calculated using the G-statistics feature in GenAlEx v6.51b2 [47] with 9999 permutations. Since the ENA-corrected FST values for null alleles computed by FreeNA were of the same magnitude as the raw values, the statistical significance (p-values) between population pairs was determined based on an exact test performed on the raw dataset in GenAlEx with 9999 permutations. We used GenAlEx ver. 6.5 [47] to estimate allele frequency by locus and population, the number of alleles (Na), the number of effective alleles (Ne), observed and expected heterozygosity values (HO and HE, respectively).
To assess potential admixture between different wildcat populations in Hungary, as well as between wildcat populations and domestic cats, we used the Bayesian clustering method implemented in version 2.3.4 of STRUCTURE [48] to determine the most likely number of genetic clusters, without predefining the populations. The Bayesian clustering method and Markov Chain Monte Carlo (MCMC) simulation were run using an admixture model and correlated allele frequencies. STRUCTURE was run with simulations of 10 runs for each value of K, with a burn-in period of 250,000 iterations and 750,000 replications for a number of genetic clusters (K) from one to seven.
The optimal number of clusters (K) was determined using the StructureSelector software (https://lmme.ac.cn/StructureSelector/, accessed on 6 October 2026) [49]. This involved calculating the second-order rate of change in the log probability of data (ΔK), as well as evaluating the four K-estimator statistics applying threshold values between 0.5 and 0.8 [50]. Following the determination of the most likely number of genetic clusters (K) using StructureSelector, the results of the 10 independent replicate runs for each K value were aligned and averaged using the CLUMPAK pipeline (https://lmme.ac.cn/StructureSelector/, accessed on 6 October 2026) [51] integrated within the StructureSelector platform. We applied the LargeKGreedy (https://lmme.ac.cn/StructureSelector/, accessed on 6 October 2026) algorithm with 2000 random input orders to compute the consensus individual membership proportions q(i). These consensus q(i) values were subsequently used to generate the final bar plots. For the K = 2 solution, field samples were classified as wildcat-associated when the membership coefficient to the wildcat cluster was q(i) ≥ 0.75, as domestic cat-associated when q(i) ≤ 0.25, and as putatively admixed for intermediate values (0.25 < q(i) < 0.75). This interval is centred on the value of q(i) = 0.5 expected for an F1 genotype and corresponds to the range conventionally applied in wildcat admixture studies [27].
To investigate the spatial genetic variability and patterns of the Hungarian wildcat population, a spatial Principal Component Analysis (sPCA) was performed [52]. For this analysis, individuals identified as domestic cats (q(i) ≥ 0.75 assignment to the domestic cat cluster) based on the STRUCTURE results were excluded from the free-living samples, leaving only putatively admixed and wildcat individuals. The analysis was carried out in the R software 4.5.0 (R Development Core Team, 2020) environment using the adegenet (version 2.1.11) and spdep (version 1.4.2) packages [52], based on data from 14 microsatellite loci. Unlike traditional Principal Component Analysis (PCA), sPCA is a multivariate method that explicitly incorporates the geographic coordinates of individuals alongside their genetic profiles. To calculate spatial autocorrelation (Moran’s I), a distance-based connection network was defined. To ensure a continuous network where each sampling point has at least one neighbour, the threshold distance was determined based on the maximum distance between nearest neighbours. During the construction of the connection network, the minimum Euclidean distance was set to 0. The maximum threshold distance was calculated dynamically, based on the maximum of the nearest-neighbour distances between individuals (k = 1 neighbour criterion), thereby ensuring spatial continuity and that every examined individual had at least one neighbour during the analysis. The sPCA identified synthetic components that maximize the product of genetic variance and spatial autocorrelation, summarizing the spatial genetic structure of the population. For the visualization of global patterns, the values of the first global principal components were represented in geographic scale along a continuous colour scale. To statistically test for Isolation by Distance (IBD), a Mantel test [53] was performed with 999 permutations using the ade4package [54] in the R 4.5.0 statistical environment. This test examined the correlation between the genetic distance matrix (Euclidean distance) and the geographic distance matrix to determine the strength (r) and statistical significance (p) of the relationship.
To evaluate the statistical power of the STR marker set for detecting hybridisation between domestic and wildcats, we simulated parental and putatively admixed genotypes using recom-sim.py [55]. As baseline data for these simulations, we selected complete empirical microsatellite profiles (i.e., featuring no missing loci) that were assigned to either parental group (wildcat-associated or domestic-cat-associated) with a membership proportion q(i) ≥ 0.75 in the initial STRUCTURE analysis. For the K = 2 solution, field samples were classified as wildcat-associated when the consensus membership coefficient to the wildcat cluster was q(i) ≥ 0.75, as domestic-cat-associated when q(i) ≤ 0.25, and as putatively admixed for intermediate values (0.25 < q(i) < 0.75). This interval is centred on the value of q(i) = 0.5 expected for an F1 genotype, and the 0.75 threshold follows established practice in European wildcat admixture studies [27]. From the free-living samples, 84 individuals were assigned to the wildcat cluster based on this assignment probability; these, along with 30 domestic cats, were used as reference populations to generate the various genotype classes. We simulated a total of 400 genotypes across the following categories: two parental classes, F1 (parental × parental), F2 (F1 × F1), first-generation backcrosses (offspring of a first-generation hybrid × parental) and second-generation backcrosses (offspring of a first-generation backcross × parental) to both parental lineages (50 simulated genotypes per class). A subsequent STRUCTURE analysis was then performed on these simulated parental and putatively admixed genotypes. Consistent with the previous run, we applied the admixture model with correlated allele frequencies, utilizing a burn-in period of 250,000 iterations followed by 750,000 Markov Chain Monte Carlo (MCMC) iterations. In this step, however, we restricted the tested number of genetic clusters to K = 1 and K = 2, performing 10 independent runs for each K value. Within each K value, the most likely of the 10 independent replicates was chosen based on the Estimated Ln Prob of Data, with the highest value representing the optimal run.
2.3.2. Mitochondrial Analysis
The ND5 sequences of the 118 field samples were checked and corrected in BioEdit version 7.7.1.0 [56]. Nucleotide sequences were aligned by ClustalW implemented in MEGA12 (version 12.0.14) [57], and manually adjusted if needed. The alignment was trimmed to 669 bp in length (positions 13,243-13,911) to be comparable with the haplotypes identified by [35], as described in their Supplementary Table S1. The number of distinct haplotypes, haplotype and nucleotide diversities, and the number of polymorphic sites were computed with DnaSP version 6.12.03 [58]. Neutrality statistics (Tajima’s D and Fu and Li’s F*) were calculated using the pegas (version 1.4) package [59] in the R environment (version 4.5.0; R Development Core Team, 2020). For the Hungarian sequences, the best nucleotide substitution model was selected by the Maximum-Likelihood analysis implemented in MEGA12 using the Bayesian Information Criterion (BIC) as suggested [57]. Relationships between the haplotypes were also calculated and visualized by constructing a haplotype network using the median-joining method available in PopART version 1.7 [60], setting ε = 10 and transversions/transition weighting of 3:1.
To investigate the hierarchical distribution of genetic variance, an analysis of molecular variance (AMOVA) was computed using the poppr (version 2.9.8) package [61] in the R environment (version 4.5.0; R Development Core Team, 2020). The analysis was run for 10,000 permutations to evaluate the significance of each parameter. For a hierarchical analysis, the populations were grouped according to regions: TD-S and TD-N into Transdanubia, NHM-Bö, NHM-Cs, NHM-M and NHM-Bü into North Hungarian Mountains and GHP into Great Hungarian Plain.
3. Results
The detailed geographic assignments, the genetic admixture proportions (q(i) values) derived from the nuclear microsatellite analysis, and the identified mitochondrial haplotypes and haplogroups of the field wildcat individuals investigated in this study are summarized in Supplementary Table S1.
3.1. Microsatellites
For both domestic cats and wildcats, all investigated loci were polymorphic. In the wildcat groups, the lowest per-locus allele count was three (locus FCA220 in the NHM-M group) and the highest was 14 (loci FCA223 and F115, both in the GHP group); averaged over all 14 loci, the mean number of alleles per group ranged from 7.07 to 10.79. The extremes occur at different loci—within the NHM-M group alone, FCA220 carried three alleles and FCA096 carried 13—so the per-locus and mean values are not contradictory. Mean observed heterozygosity (HO) was 0.73–0.82 (in the North Hungarian Mountains—Börzsöny and Transdanubian-South groups), while mean expected heterozygosity (HE) was 0.76–0.82 (in the North Hungarian Mountains—Bükk and Great Hungarian Plain groups). Observed heterozygosity (HO) was higher than expected heterozygosity (HE) only in Southern Transdanubia (TD-S) and the Bükk group (NHM-Bü), whereas HE was higher than HO across all remaining wildcat groups. The effective number of alleles ranged from 4.56 to 5.88 (Supplementary Table S3). Departures between HO and HE at individual loci were bidirectional—54% of the locus × group combinations showed HO < HE and 46% showed HO > HE—and averaged over the 14 loci the two statistics were closely similar in every group, with a mean heterozygote deficit across the wildcat groups of 3.5%. A repeatable deficit was confined to a small number of loci, notably FCA220, consistent with the null alleles expected when markers developed for the domestic cat are applied to the wildcat; these were screened with MICRO-CHECKER and the population comparisons were repeated with the ENA correction implemented in FreeNA, which gave FST values of the same magnitude as the uncorrected estimates. Treating spatially structured groups as single populations (a Wahlund effect, consistent with the significant differentiation among populations within regions, ΦSC = 0.083, p = 0.0226) will also contribute to the small positive values observed.
The STRUCTURE program showed the highest average log-likelihood values for two genetic clusters, and the second-order rate of change in the log-likelihood values was also highest for two genetic clusters, K = 2 (Figure 2A). In this case, domestic cat references clustered into one group, while the various wildcat groups clustered into another. Within the geographically defined wildcat groups, 18 individuals were genetically identified to the domestic cat cluster, and 16 individuals were identified as putatively admixed. The four supervised statistics of the Puechmaille method indicated the presence of three clusters across all thresholds, K = 3 (Figure 2B). In this case, domestic cats and wildcats formed distinct genetic clusters; however, the internal genetic structure observed within the wildcat group (K = 3) did not correlate with geographical location. Because these internal wildcat genetic lineages do not form geographically segregated populations across the study area, all subsequent population genetic analyses were carried out under the K = 2 model, assuming that wildcats represent a single polymorphic population. Within the geographically defined wildcat groups, domestic cat (n = 18) and putatively admixed individuals (n = 10) were also identified under this model (Figure 2B).
Based on the results of the Structure clustering analysis, the NHM-Bü sub-region exhibited the lowest number of putatively admixed and domestic cat individuals among the investigated geographic groups (n = 0 in both categories). In contrast, the highest number of putatively admixed individuals was identified in the GHP population (n = 7; 20%), which also presented the highest occurrence of domestic cats (n = 6; 17.14%). Examining the further areas, out of the 9 samples from the TD-N region, 2 (22.22%) were genetically assigned as wildcats, 4 (44.44%) as putatively admixed individuals, and 3 (33.33%) as domestic cats. The vast majority of the 14 animals originating from the TD-S region were wildcats (n = 12; 85.71%), alongside which only 1 putatively admixed individual (7.14%) and 1 domestic cat (7.14%) were identified. Among the 13 analyzed specimens from the NHM-Bö group, 9 (69.23%) were classified as wildcats, while 2 individuals each (15.38%) were assigned as putatively admixed individuals and domestic cats. Within the 14 samples from the NHM-Cs group, 11 (78.57%) resulted as wildcats, 1 as a putatively admixed individual (7.14%), and 2 as domestic cats (14.29%). Finally, out of the 15 individuals in the NHM-M region, 10 (66.67%) represented the wildcat lineage, while an additional 1 putatively admixed individual (6.67%) and 4 domestic cats (26.67%) were detected (Figure 3).
Based on the results of the spatial Principal Component Analysis (sPCA), the genetic variability of the Hungarian wildcat population is not distributed randomly across space. A spatial genetic structure can be visually identified in Figure 4. Spatially close individuals (indicated by points with similar colour shades) exhibit higher genetic similarity, whereas distinct genetic differentiation is evident among populations located in more distant regions of the country. This visually identified spatial structure was statistically corroborated by the Mantel-test results. A significant positive correlation was detected between genetic and geographic distances (r = 0.118, p = 0.005). This association is statistically significant but weak, and is consistent with a pattern of isolation by distance (IBD) within the sampled population rather than demonstrating it conclusively.
Genetic distances among the investigated populations were estimated using pairwise FST values corrected for the presence of null alleles (FreeNA) (Table 1). The results revealed that the domestic cat population is strongly and highly significantly differentiated from all other investigated groups (FST = 0.034–0.139; p < 0.01 and p < 0.001). The highest degree of differentiation within the entire dataset was observed between the domestic cat and the NHM-Bü population (FST = 0.139, p < 0.001), most likely because no domestic cat or putatively admixed individuals were identified within this group. Genetic structuring among the non-domestic groups was more moderate but well-defined. The NHM-Bü population exhibited statistically significant, moderate differentiation from the other groups (FST = 0.026–0.062; p < 0.05–0.001). In contrast, among certain geographic or categorized population pairs—particularly among the NHM-Bö, NHM-Cs, and NHM-M groups—pairwise FST values were not significantly different from zero (FST ≤ 0.017; p > 0.05; negative estimates were treated as zero). No significant differentiation was therefore detected among these groups, which is consistent with gene flow between them; given the small and unequal group sizes, however, this result should not be read as positive evidence for the absence of barriers to gene flow.
The off-diagonal elements represent genetic distances between the specified population groups. Asterisks indicate the level of statistical significance: *** p < 0.001, ** p < 0.01, * p < 0.05, while “ns” represent non-significant differentiation (p > 0.05). GHP (Great Hungarian Plain), TD-N (Transdanubia North), TD-S (Transdanubia South), NHM-Bö (North Hungarian Mountains—Börzsöny), NHM-Cs (North Hungarian Mountains—Cserhát), NHM-M (North Hungarian Mountains—Mátra), and NHM-Bü (North Hungarian Mountains—Bükk).
Using the previously applied assignment probability of q(i) ≥ 0.75, 100% of the generated wildcat genotypes assigned to one cluster, and 100% of the domestic cat genotypes assigned to the other cluster. Among the generated F1s and F2s, 100% and 90% of the genotypes, respectively, were identified as putatively admixed individuals; none of the F1 genotypes was assigned to a parental cluster, whereas the remaining 10% of the F2 genotypes were assigned to the domestic cat cluster. Of the generated second backcrosses, 76% and 84% assigned to the corresponding parental cluster, whereas 24% and 16% failed to assign. 100% of the third backcrosses assigned to the appropriate parental group, with no generated genotypes remaining unassigned. Detailed information regarding the Structure q(i) values for the generated genotypes is provided in Supplementary Tables S4 and S5.
3.2. Phylogenetic Analysis
By analyzing the 669 bp sequence portion, we identified eight haplotypes among the 118 samples. The aligned sequences contained 16 polymorphic sites, corresponding to 16 mutations; under an infinite-sites model without recombination this number of segregating sites can support up to 17 distinct haplotypes, of which eight were observed here. The total haplotype diversity (Hd) of the analyzed dataset was 0.683 ± 0.036, and the nucleotide diversity (π) was 0.00704 ± 0.00062. The average number of nucleotide differences (k) among individuals was 4.712. Demographic statistics for the overall dataset did not show significant deviations from the mutation-drift equilibrium, yielding a Tajima’s D value of 1.5719 (p = 0.1160) and a Fu and Li’s F* value of 0.6815 (p = 0.4617) (Table 2).
The eight identified haplotypes belonged to four main haplogroups: sub-clade Ia (n = 60; 50.85%) was represented by the dominant haplotype dw4, followed by haplogroup III (n = 25; 21.19%, driven mainly by haplotype dw1, n = 22; 18.64%), haplogroup V (n = 15; 12.71%; w2), haplogroup IV (n = 12; 10.17%; w1), and sub-clade Ib (n = 6; 5.08%; dw3), while haplotype dw5 was found in only a single sample (0.85%). Two new haplotypes (Hap1 and Hap7) were identified, both belonging to haplogroup III and each represented by a single sample (0.85%).
When the dataset was partitioned into the seven geographic sub-regions, the number of identified haplotypes ranged from three (TD-S, NHM-Cs, and NHM-Bü) to five (GHP) (Figure 5). Haplotype dw4 predominated across the NHM populations (NHM-Bö, NHM-Cs, NHM-M, and NHM-Bü), whereas dw1 was most common in the GHP group, and both dw1 and dw4 were equally prevalent in TD-N. In contrast, the TD-S population was characterized by a predominance of haplotype w1. Genetic variability was highest in the TD-N and GHP sub-regions, both in terms of haplotype and nucleotide diversity, while the lowest diversity values were observed in the NHM-Cs and NHM-Bö groups (Table 2). Demographic statistics (Tajima’s D and Fu and Li’s F*) did not reveal significant deviations from equilibrium in any sub-region (p > 0.05), although the NHM-Bö group showed a negative trend approaching significance (Tajima’s D = −1.7203, p = 0.0854; Fu and Li’s F* = −1.9014, p = 0.0699).
3.3. Phylogeographic Analysis of Hungarian Wildcat and Putative Admixed Individuals
The reconstructed MJ network, which includes the eight haplotypes (dw1; dw3; dw4; dw5; w1; w2; Hap1; Hap7) identified from the analyzed dataset, revealed a distinct genetic structure. The network topology is structured around two highly frequent and dominant haplotypes: dw4 and dw1. These two central variants are separated by a total of three mutational steps through an intermediate, lower-frequency haplotype (dw3), which is separated by one mutation from dw4 and two mutations from dw1. Node dw1 forms the centre of a star-like radiation: low-frequency, likely private or rare haplotypes (Hap1, Hap7, dw5) diverge from this central variant by a single mutational step. Conversely, w1 and w2 represent a genetically more distant, divergent branch. w1 is separated from the central dw1 variant by seven mutational steps, followed by three additional mutations leading to haplotype w2. The colour-coded geographic repartition highlights that the most frequent haplotypes (particularly dw4 and dw1) are not strictly isolated to a single geographic region but are shared among several sub-regions (e.g., the Great Hungarian Plain/GHP and various North Hungarian Mountains/NHM groups). This distribution pattern supports the presence of genetic connectivity and historical gene flow among the populations (Figure 6).
To investigate the hierarchical distribution of genetic variance, an analysis of molecular variance (AMOVA) was performed, revealing that the great majority of the variation was partitioned within samples (σ2 = 3.899; 77.88% of the total variance; ΦST = 0.221; p = 0.0001). Differentiation among populations within regions accounted for a small but statistically significant share of the variance (σ2 = 0.354; 7.08%; ΦSC = 0.083; p = 0.0226). In contrast, differentiation among the three geographic macro-regions explained 15.04% of the variance and was not statistically significant (σ2 = 0.753; ΦCT = 0.150; p = 0.1545) (Table 3).
4. Discussion
The microsatellite loci used in the genetic analysis of the Hungarian wildcat population—the majority of which were originally developed for the domestic cat [43]—performed adequately for estimating population genetic metrics. The average number of alleles (NA: 7.07–10.79) and the expected heterozygosity (HE: 0.76–0.82) measured in the Hungarian samples are high. This allelic diversity exceeds the values reported for several other European populations, such as those in Germany (HE: 0.65) and Italy (HE: 0.69) [62,63]. The Hungarian results are most comparable to values observed in France (NA: 7.0; HE: 0.73) and Croatia (NA: 7.2; HE: 0.72) [64,65]. It is crucial to note, however, that the seemingly higher diversity indices observed in wildcat samples collected from natural habitats often do not reflect population “health,” but rather suggest the presence of historical or recent introgression (gene flow) from domestic cats, as demonstrated by [33] in the case of Swiss populations.
Although STRUCTURE analysis (at K = 3) identified clusters within the free-living wildcat samples, these did not exhibit any distinct spatial patterns. These results suggest that genetic differentiation between certain regions (e.g., Börzsöny, Cserhát, Mátra) is limited, which suggests weak regional differentiation, substantial gene flow, and continuous admixture. The East European—and specifically the Hungarian—wildcat population is genetically highly admixed and integrated; thus, it cannot be differentiated into distinct biogeographic units, as has been documented in populations such as those in Scotland. This may result from long-standing hybridisation and introgression [4]. While clear genetic discontinuities caused by mountain ranges or major highways have been observed in other European populations [62,66], such significant natural or anthropogenic barriers to gene flow are absent among the study areas in Hungary.
In contrast, sPCA revealed high-resolution structure, as this method incorporates the geographic coordinates of the samples. Following the removal of identified domestic cat individuals, the Mantel-test suggest a significant but weak positive correlation between genetic and geographic distances (r = 0.118, p = 0.005). This suggests that isolation by distance is present. The Hungarian population is characterized not by strong genetic breaks, but by gradual, continuous transitions (clines), which may suggest well-connected habitats and high mobility of the species [66].
Based on these results, it can be stated that hybridisation between wildcats and domestic cats occurs in Hungarian natural populations. In our current study, we detected an intermediate assignment rate of 13.56% based on the STRUCTURE analysis, while free-ranging domestic cats identified in the field accounted for an additional 15.25% of the samples (bringing the combined non-wildcat proportion to 28.81%). As a methodological limitation, it should be noted that the parental wildcat reference individuals for the simulations were likewise selected based on a STRUCTURE q(i) > 0.75 threshold rather than being independently validated pure wildcat references derived from genome-wide genomic methods. This is consistent with findings in various other European populations; however, the extent of hybridisation appears to vary significantly across regions. While hybridisation rates in Germany are relatively low, estimated at approximately 3% [27,32], higher levels indicating moderate hybridisation have been documented in other countries, such as Italy (3–11%), Switzerland (12.5–15%), and the Iberian Peninsula (14–21%) [11,18,32,37,62,67]. Refs. [64,65] reported hybridisation rates of approximately 29% in France, with figures ranging from 13% to 52% recorded across Balkan populations. The most extreme rates of hybridisation have been observed in Scotland, where values can reach 80–90% [12,24]. In Hungary, preliminary studies conducted in the mid-2000s estimated hybridisation rates between 25% and 31% using genetic methods [4,30].
During the 18th and 19th centuries, the European wildcat was subjected to intensive hunting, leading to drastic population declines, local extinction, and fragmentation [68]. When population densities reach critical levels, the scarcity of conspecifics increases the likelihood of mating with the ubiquitous domestic cat [4]. Deforestation and habitat fragmentation have forced wildcats into small, isolated habitat [68]. This contributes to a decline in genetic diversity and heightens the risk of hybridisation, particularly at habitat borders [30]. The mosaic of agricultural areas and forest patches (agro-forestry ecosystems) facilitates the interaction of the domestic cat and wildcat. While continuous forests provide refuge for wildcats, settlements and agricultural lands serve as centres for domestic cats and hybrids; within these contact zones, hybridisation becomes more frequent [17,24]. The extensive hybridisation observed in Hungarian wildcats was historically predestined by habitat loss and population declines resulting from persecution. This demographic crisis, compounded by the mosaic landscape structure characteristic of Hungary, enables increased contact between wild and domestic cats along forest edges. Together with our current genetic findings, these factors substantiate the ancient and persistent nature of introgression, which basically threatens the genetic integrity of the local population [4,30].
In addition to the loss of genetic purity, this ongoing interbreeding can have serious ecological consequences. As a mesocarnivore, the wildcat provides essential top-down control of small mammal populations within forest ecosystems [22]. Extensive hybridisation may alter phenotypic traits, hunting behaviour, and habitat use, potentially shifting the ecological niche of admixed individuals closer to that of feral cats, thereby disrupting their native ecosystem function [14]. Therefore, strict and targeted conservation management is urgently required. Since our spatial analyses (sPCA) revealed gradual clines rather than sharp genetic breaks, management strategies must focus on the agro-forestry contact zones. To mitigate further introgression, conservation strategies in various European regions often focus on establishing buffer zones around core forest habitats [6,32]. Within these contact zones, wildlife management practices frequently involve the regulation or removal of free-roaming and feral domestic cats [32]. However, since morphological identification has proven unreliable in highly admixed populations [31], such removal efforts in several countries are increasingly supported by rapid genetic screening to avoid the accidental culling of pure wildcats [33]. Furthermore, international conservation frameworks highlight that public awareness campaigns targeting rural cat owners regarding responsible pet ownership (e.g., sterilization and indoor housing) serve as essential components of comprehensive management plans aimed at reducing the constant influx of domestic genes into the wildcat gene pool [12,27].
The mitochondrial DNA analysis of the Hungarian wildcat population revealed the presence of eight distinct haplotypes, which form a structured median-joining network centred around two dominant variants, dw4 and dw1. Evaluated within the broader European phylogeographic context, these haplotypes are clearly split into two deeply divergent evolutionary lineages: Lineage D, which encompasses the domestic cat and African wildcat origins, and Lineage W, which represents wildcat-assigned ancestry [10,35,36]. The widespread distribution of shared (dw) haplotypes across various Hungarian sub-regions initially appears to confirm extensive, country-wide introgression of domestic cats. However, a detailed phylogenetic assessment reveals a much more complex evolutionary history. One of the most notable phylogeographic patterns in the Hungarian populations the overwhelming dominance of the dw4 haplotype in the area of the North Hungarian Mountains (NHM-Bö, NHM-Cs, NHM-M, and NHM-Bü). This variant belongs to sub-clade Ia of lineage D, which, according to evolutionary models, diverged approximately 37,700 years ago [35]. Since this divergence long pre-dates currently known archeological and genetic evidence for cat domestication [10,69,70], the high frequency of dw4 in northern habitats does not necessarily suggest recent introgression with modern domestic cats. Rather, one possible explanation is a legacy of late Pleistocene, historical lineage retention and ancient introgression, when European wildcats migrating toward the Levant and Anatolia interbred with local Felis lybica populations [36,71,72]. This ancient introgression may reflect the mito-nuclear discordance observed in the regions of the North Hungarian Mountains, where individuals appear genetically pure based on nuclear microsatellites yet carry a shared mitochondrial pattern [26]. In contrast, the dw1 haplotype, predominantly found in the Great Hungarian Plain (GHP), reveals an entirely different evolutionary story. Based on demographic simulations, haplotypes such as dw1 and dw3 differentiated recently, about 1000 years ago, and may suggest recent human-mediated introgression originating from Felis catus [35]. This recent introgression is in perfect concordance with the high proportion of admixed individuals detected by nuclear STR analysis in the Great Hungarian Plain, reflecting intense anthropogenic pressure and frequent interactions with domestic cats in highly fragmented agro-forestry mosaic landscapes [4,6]. The ancestral evolutionary heritage of the European wildcat (Lineage W) in the Hungarian population is primarily preserved by the w1 and w2 haplotypes. Notably, the w1 haplotype is particularly frequent in Southern Transdanubia (TD-S), a region that also exhibits a high proportion of wildcat-assigned individuals based on nuclear results. According to the continental phylogeographic framework, the w1 variant belongs to clade IV, corresponding to a South-Eastern European spatial cluster [35]. The dominance of this wildcat-assigned lineage in southern Hungary likely reflects a post-glacial recolonization wave originating from Balkan or extra-Mediterranean glacial refugia following Pleistocene climatic oscillations [4,29,73]. The complex yet continuous nature of these genetic patterns is further supported by the analysis of molecular variance (AMOVA). The hierarchical partitioning of genetic variance showed that the vast majority of variation (77.88%) is found within individual samples, whereas genetic differentiation among the main geographic regions (15.04%) was not statistically significant. The absence of sharp, large-scale regional isolation correlates with the topology of the median-joining network, which clearly illustrates that dominant haplotypes are heavily shared across geographically distant sub-regions. Overall, these results suggest that the Hungarian wildcat population has not been characterized by long-term geographic isolation; its genetic structure has been continuously shaped by the species’ high mobility, historical recolonization waves, ancient inter-taxon gene flow, and varying local intensities of recent anthropogenic hybridisation [35,66].
5. Conclusions
Our study provides a comprehensive overview of the genetic diversity, spatial structure, and ecological challenges posed by hybridisation in the Hungarian European wildcat population. The patterns compatible with historical lineage retention and ancient introgression observed in the North Hungarian Mountains (presumably originating from past admixture with the Felis lybica lineage) and the recent hybridisation with domestic cats detected in the Great Hungarian Plain well reflect the complex evolutionary and demographic history of the population. Spatial analyses (sPCA) did not reveal sharp genetic barriers among the sampled Hungarian habitats; the population is characterized by continuous spatial transitions and isolation by distance, consistent with spatial gene flow across the landscape. Nevertheless, extensive hybridisation severely threatens the genetic integrity of wildcats and their ecological role in regulating small mammal populations within forest ecosystems. Identification based on morphological traits is no longer reliable due to the high degree of admixture. Therefore, future conservation strategies must focus on establishing buffer zones in critical agro-forestry contact zones, supporting the regulation of feral domestic cats through rapid and reliable genetic screening methods, reconsidering the current legal regulation of naturally occurring putatively admixed individuals, and launching public awareness campaigns on responsible pet ownership (e.g., spaying/neutering, indoor housing), which are essential to halt the continuous influx of domestic cat genes. The combined analysis of microsatellite loci and the mitochondrial ND5 gene suggests that the sampled population retains high genetic diversity at the markers examined, and is compatible with the presence of both historical and recent introgression. Given the sample size and the resolution of the marker panel, these patterns should be regarded as indicative rather than conclusive. In particular, diversity estimates for sub-regions with small sample size (e.g., TD-N, n = 9) are subject to substantial sampling error and should be interpreted with caution rather than used as the basis for direct regional comparisons.