Disease-driven mass mortality event leads to widespread extirpation and variable recovery potential of a marine predator across the eastern Pacific

AbstractThe prevalence of disease-driven mass mortality events is increasing, but our understanding of spatial variation in their magnitude, timing, and triggers are often poorly resolved. Here, we use a novel range-wide dataset comprised of 48,810 surveys to quantify how Sea Star Wasting Disease affected Pycnopodia helianthoides, the sunflower sea star, across its range from Baja California, Mexico to the Aleutian Islands, USA. We found that the outbreak occurred more rapidly, killed a greater percentage of the population, and left fewer survivors in the southern half of the species’ range. Pycnopodia now appears to be functionally extinct (> 99.2% declines) from Baja California, Mexico to Cape Flattery, Washington, USA and exhibited severe declines (> 87.8%) from the Salish Sea to the Gulf of Alaska. The importance of temperature in predicting Pycnopodia distribution rose 450% after the outbreak, suggesting these latitudinal gradients may stem from an interaction between disease severity and warmer waters. We found no evidence of population recovery in the years since the outbreak. Natural recovery in the southern half of the range is unlikely over the short-term and assisted recovery will likely be required for recovery in the southern half of the range on ecologically-relevant time scales., MethodsThirty research groups from Canada, the United States, Mexico, including First Nations, shared 34 datasets containing field surveys of Pycnopodia (Table S1). The data included 48,810 surveys from 1967 to 2020 derived from trawls, remotely operated vehicles, SCUBA dives, and intertidal surveys. We compiled survey data into a standardized format that included at minimum the coordinates, date, depth, area surveyed, and occurrence of Pycnopodia for each survey. When datasets contained more than one survey at a site in the same day (e.g. multiple transects), we divided the total Pycnopodia count in all surveys by the total survey area and averaged the latitude, longitude, and depth as necessary. Using breaks in data coverage, political boundaries, and biogeographic breaks we assigned each survey to one of twelve regions: Aleutian Islands, west Gulf of Alaska (GOA), east Gulf of Alaska, southeast Alaska, British Columbia (excluding the Salish Sea), Salish Sea (including the Puget Sound), Washington outer coast (excluding the Puget Sound), Oregon, northern California, central California, southern California, and the Pacific coast of Baja California (Fig. S1; see Supplementary Material)., Usage notesDocumentation, data, and code accompanying Hamilton et al., 2021 Pycnopodia Rangewide Assessment paper. Data MasterPycno_ToShare: Dec_lat = latitude in decimal degrees. Numeric. Dec_lon = longitude in decimal degrees. Numeric. Depth = depth in meters. Numeric. Pres_abs = presence or absence of Pycnopodia on that survey. Binary. Presence = 1, absence = 0 Density_m2 = density in meters squared if available for that set of surveys. Numeric. NA = no density data available for that survey. Source = shorthand name of the group that shared the data with us and the type of data (e.g. trawl, dive). To get further info on who that dataset, group, and group contact, see Table S1. Character. Note: When datasets contained more than one survey at a site in the same day (e.g. multiple transects), we divided the total Pycnopodia count in all surveys by the total survey area and averaged the latitude, longitude, and depth as necessary in order to minimize the impacts of pseudoreplication on the dataset. Used in MaxentSWD_Final and Density-Inc_Models_Figs_Tables_ToShare. CrashEventsForRPlot: Crash Dates were determined trends in Pycnopodia occurrence (site-level presence or absence) to estimate ‘crash date’, defined as the date when the occurrence rate of Pycnopodia in a region decreased by 75% from pre-outbreak levels. Used in OutbreakTimelineFigs_ToShare.R EpidemicPhases: See manuscript methods for information on how the column ‘EpidemicPhases’ was created. “Start-End” specifies whether that date was the start or the end of that epidemic phase for that region. Used in OutbreakTimelineFigs_ToShare.R Incidence_2012-2019: Columns G-J were calculated by fitting a logistic regression model to the occurrence of Pycnopodia over time for each region. We fit a logistic regression model to the occurrence of Pycnopodia from 1/1/2012 to 12/31/2019 to model the shape of the population decline for each region (Fig. 1a). From these models, we 1) estimated regional Pycnopodia occurrence rates on 1/1/2012 and 12/31/2019, 2) calculated the predicted occurrence value corresponding to a 75% decline in starting versus ending occurrence in each region, and 3) solved the inverse logistic equations for the date at which this occurrence value was predicted. All other columns are identifying information derived from MasterPycno_ToShare. Used in OutbreakTimelineFigs_ToShare.R MasterPycno_021821_SpatialJoin: Used to make Fig 5 for the remnant population analysis. Grid IDs were created in ArcGIS using a spatial join. Used in RemPop_Final.R Maxent_SWD OccurrenceData, BiasFiles, Maxent_EnvironmentalFiles folders all contain data that is used to run Maxent Analysis in Maxent terminal (available here: https://biodiversityinformatics.amnh.org/open_source/maxent/) Maxent_Output contains the output files generated by the Maxent terminal. Used in Maxent_SWD_Final to make Figures 3cd and 4. Code Density_Inc_Models_Figs_Tables_Toshare: Makes Fig 2, runs models on density and incidence changes across the phases of the pandemic ContinuedDeclines_Models_Figs_Tables_Toshare: Makes Fig S3, runs models on density and incidence changes for the years 2017-2020 (generally after the outbreak had subsided in many regions). OutbreakTimelineFigs_ToShare: Makes Fig1 Maxent_SWD_Final: Prepares Pre and Post-outbreak Pycnopodia incidence and bias files for Maxent program. Maxent analysis must then be run using these files and the Maxent program. For this analysis, we were not able to use Dismo or other R packages that run Maxent because these packages do not give you the ability the adjust the parameter Tau. Maxent may be downloaded at https://biodiversityinformatics.amnh.org/open_source/maxent/. Using the Maxent output files (included in the directory), this code then creates Figure 4 from the manuscript. Figure 3 were made using density data and Maxent output visualized in ArcGIS. RemPop_Final.R: Creates Figure 5, runs remnant population analysis. Uses MasterPycno_021821_SpatialJoin.摘要 疾病驱动的大规模死亡事件发生率持续上升，但我们对其规模、发生时间与触发因素的空间差异往往尚不清楚。本研究利用一项包含48810次调查的全新全分布区数据集，量化了海星消瘦病(Sea Star Wasting Disease)对向日葵海星（Pycnopodia helianthoides）的影响，其分布范围覆盖从墨西哥下加利福尼亚州到美国阿留申群岛的整个区域。研究发现，疫情在该物种分布区南部的传播速度更快、种群死亡率更高、幸存个体更少。目前，从墨西哥下加利福尼亚州到美国华盛顿州弗拉特利角的区域内，向日葵海星已出现功能性灭绝（种群数量下降超99.2%）；从萨利希海至阿拉斯加湾的区域内，其种群也出现了严重下降（下降幅度超87.8%）。疫情暴发后，温度在预测向日葵海星分布的重要性提升了450%，这表明上述纬度梯度可能源自病害严重程度与暖水之间的交互作用。研究未发现疫情暴发后种群出现恢复的迹象。分布区南部的种群在短期内难以自然恢复，在生态相关的时间尺度上，南部种群的恢复可能需要人为辅助恢复措施。 方法 来自加拿大、美国、墨西哥（含原住民群体）的30个研究团队共享了34份包含向日葵海星野外调查数据的数据集（表S1）。该数据集涵盖1967年至2020年间的48810次调查，数据来源包括拖网、遥控水下机器人（remotely operated vehicle, ROV）、水肺潜水（SCUBA）以及潮间带调查。我们将所有调查数据整合为标准化格式，每份记录至少包含坐标、调查日期、水深、调查区域面积以及向日葵海星的出现情况。当数据集包含同一地点单日多次调查（如多条样带）时，我们将所有调查中的向日葵海星总数量除以总调查面积，并在必要时对纬度、经度和水深进行平均处理。基于数据覆盖断点、行政边界和生物地理分区断点，我们将每份调查划分为12个区域之一：阿留申群岛、阿拉斯加湾西部（west Gulf of Alaska, GOA）、阿拉斯加湾东部、阿拉斯加东南部、不列颠哥伦比亚省（不含萨利希海）、萨利希海（含普吉特湾）、华盛顿州外海海岸（不含普吉特湾）、俄勒冈州、北加利福尼亚州、中加利福尼亚州、南加利福尼亚州以及下加利福尼亚州太平洋沿岸（图S1；详见补充材料）。 使用说明 本数据集配套Hamilton等人2021年发表的《向日葵海星全分布区评估》论文的文档、数据与代码。 数据集MasterPycno_ToShare： Dec_lat：十进制纬度，数值型。 Dec_lon：十进制经度，数值型。 Depth：水深（单位：米），数值型。 Pres_abs：本次调查中向日葵海星的出现与否，二值型，出现为1，未出现为0。 Density_m2：单位面积种群密度（单位：平方米），若该批次调查有密度数据则为数值型，无密度数据则为NA。 Source：共享数据的团队简称与数据类型（如拖网、潜水），字符型。如需获取该数据集所属团队及联系方式的详细信息，请参见表S1。 注：当数据集包含同一地点单日多次调查（如多条样带）时，我们将所有调查中的向日葵海星总数量除以总调查面积，并在必要时对纬度、经度和水深进行平均处理，以尽可能降低伪重复对数据集的影响。该数据集用于MaxentSWD_Final与Density-Inc_Models_Figs_Tables_ToShare。 CrashEventsForRPlot：通过分析向日葵海星出现情况（样点水平的出现/未出现数据）的趋势来确定“崩溃日期”，即某区域向日葵海星出现率较暴发前下降75%的日期。该数据集用于OutbreakTimelineFigs_ToShare.R。 EpidemicPhases：关于“EpidemicPhases”列的构建方法，请参见论文方法部分。"Start-End"用于标注该日期属于该区域疫情阶段的起始还是结束。该数据集用于OutbreakTimelineFigs_ToShare.R。 Incidence_2012-2019：G-J列通过对每个区域的向日葵海星出现情况随时间的变化拟合逻辑回归模型计算得到。我们对2012年1月1日至2019年12月31日的向日葵海星出现数据拟合逻辑回归模型，以模拟各区域的种群下降趋势（图1a）。基于这些模型，我们1）估算了2012年1月1日和2019年12月31日各区域的向日葵海星出现率；2）计算了对应各区域初始出现率较最终出现率下降75%的预测出现值；3）通过求解逆逻辑回归方程得到该预测出现值对应的日期。其余列均为源自MasterPycno_ToShare的标识信息。该数据集用于OutbreakTimelineFigs_ToShare.R。 MasterPycno_021821_SpatialJoin：用于制作残存种群分析的图5。网格ID通过ArcGIS的空间连接工具生成。该数据集用于RemPop_Final.R。 Maxent_SWD文件夹下的OccurrenceData、BiasFiles、Maxent_EnvironmentalFiles均包含用于在Maxent终端中运行Maxent分析的数据（Maxent终端下载地址：https://biodiversityinformatics.amnh.org/open_source/maxent/）。Maxent_Output文件夹包含Maxent终端生成的输出文件。该数据集用于Maxent_SWD_Final以制作图3c、d和图4。 代码文件Density_Inc_Models_Figs_Tables_Toshare：用于制作图2，运行针对疫情各阶段种群密度与出现率变化的模型。 ContinuedDeclines_Models_Figs_Tables_Toshare：用于制作图S3，运行针对2017-2020年（多数区域疫情已平息后）的种群密度与出现率变化的模型。 OutbreakTimelineFigs_ToShare：用于制作图1。 Maxent_SWD_Final：为Maxent程序准备疫情暴发前后的向日葵海星出现率与偏差文件。需使用该文件与Maxent程序开展分析。本次分析无法使用Dismo或其他运行Maxent的R包，因为这些包无法调整参数Tau。Maxent程序可从https://biodiversityinformatics.amnh.org/open_source/maxent/下载。该代码会利用目录中包含的Maxent输出文件生成论文中的图4。图3通过密度数据与Maxent输出文件在ArcGIS中可视化制作。 RemPop_Final.R：用于生成图5，开展残存种群分析，使用MasterPycno_021821_SpatialJoin数据集。