我叫大卫·伯纳姆,这次我想探讨的是,如何通过模拟产品随时间推移发生的降解过程来估算其保质期。一种方法是将产品放入温控箱中。然而,这种方法可能需要数月甚至数年的时间才能积累足够的数据来构建模型。但如果我们推出的新产品本身就不稳定呢?我们需要尽快获得这方面的信息。
加速稳定性研究的作用就是在短短几周内为我们提供这些数据。我们通过对产品施加压力来实现这一点,例如,我们可以将其储存在高温环境中。我开发了一款 JSL 应用程序,旨在对这些加速稳定性研究产生的数据进行建模。这项工作是 Pega Analytics 与几家制药公司合作完成的。Pega Analytics 是我 15 年前创立的公司,主要提供 JMP 培训,但也提供 JSL 服务。今天,我想在演讲中探讨一下药物片剂降解建模的过程。
这种降解现象的发生是因为这些药片并非完全稳定。其中含有会随时间演变的杂质。我可以用动力学模型来科学地描述这种演变过程。杂质的演变是由药片内部和环境中的热能和水分驱动的。让我们来看一下这类研究的数据。研究中使用的是高温(即储存温度)和高相对湿度。由于这种加速作用,因为我们在短短两周内对产品进行了高强度处理,所以可以测量到杂质的含量,并且这些数值的变化足以让我建立一个模型。
让我给您展示一下模型的样子。这个模型包含三个部分:温度和相对湿度,这些是储存条件;时间部分,代表杂质随时间的演变。一旦这些数据足以构建模型,我就可以做两件事。首先,我想降低储存条件,使其恢复到正常的无应力状态。然后,我想进行外推,尝试预测该产品未来的变化趋势。我调整一下Y轴的刻度,以便您更清楚地看到。在这个例子中,您可以看到时间曲线非常平缓,说明产品非常稳定,而且低于紫色线,也就是该产品1%的规格阈值。这个动力学模型看起来一切正常。
我们还可以构建另一种模型,即包装模型。让我们解释一下为什么要进行包装。假设我有一个不稳定的产品。为了说明这一点,我将温度调回 50 度。现在您可以看到,大约 50 到 60 天后,这种杂质含量超过了规格限值,因此我们的产品不再合格。想象一下,我知道在这个例子中我将温度调到了 50 度,但假设这是正常的储存条件。那么我们就需要考虑可以采取什么措施。我们可以引入包装。包装的作用是降低湿度。
随着湿度降低,你会看到时间曲线变得平缓。事实上,随着湿度降低,温度敏感性也随之降低。这里我设定的相对湿度是20%,这是我在受控环境下进行包装时可以达到的。我把药片放入瓶中,在20%的湿度下密封瓶子。其实我不需要包装模型就能做到这一点,我只需要把药片放进瓶子里就行了。但随着时间的推移,瓶内的湿度会逐渐升高。这真是太神奇了。你可能觉得瓶子很坚固,但实际上水分会从瓶壁渗透进来,我们可以对这种渗透过程进行建模。这真是一个非常有趣的物理现象。
我们可以通过添加干燥剂来减缓干燥过程,并且可以模拟干燥剂的作用。利用这些模型,我们可以根据干燥剂的用量进行假设分析。这就是我们对包装模型感兴趣的原因。虽然它非常小众,但应用范围很广,也很有趣,因为它涉及物理和数学。为了让更多人能够理解,我将重点讲解动力学模型。即使我们研究的是包装,如果我们不了解湿度对干燥过程的影响,包装模型也毫无用处。我将重点讲解这些动力学模型。我还想探讨一下构建这些模型时遇到的挑战。不是模型构建本身,而是模型选择。
我认为模型选择方面的挑战在加速稳定性研究中非常独特。JSL 应用程序将帮助我解决这个问题。但在介绍应用程序本身之前,我想先让大家了解一下我们如何使用 JMP 的原生功能构建这些模型的一些基本原理。我将向大家介绍我认为针对这些数据可以构建的最简单的模型。这里有一个模型公式。这个公式基于阿伦尼乌斯方程。实际上,它包含一个温度的指数项、一个时间项以及一个相对湿度项。我们扩展了阿伦尼乌斯方程,使其包含水分。现在我可以向大家展示这个公式的性能。
如果我取消隐藏这一列,公式就会预测我的杂质含量,但结果并不理想。它预测的是 10%,而实际值是 2.5%。这是怎么回事呢?嗯,如果我再看一下这个公式,您应该知道,公式是基于不同列变量的表达式。但您可能不熟悉,也可能了解,公式还可以添加参数。我可以添加一个参数,给它命名并赋值。现在,这个公式有三个参数,每个参数的值都是随意分配的。这不是我想要的。我想要的是为这些参数选择最佳值。
我通过将模型拟合到我的数据来实现这一点。这个过程是在 JMP 的非线性平台中完成的。接下来,我将更详细地了解一下非线性平台。首先,我将模型分配给我的预测公式,并将杂质分配给我的响应。这会弹出一个窗口,上面写着“点击开始”。在点击“开始”按钮之前,让我先解释一下它的作用。我们有一些参数,并且我已经为它们赋予了初始值。我想要做的是优化这些值。当我点击“开始”按钮时,它会应用高斯-牛顿算法,该算法旨在选择使误差平方和最小的值。这与最小二乘回归的概念类似。如果您愿意,您可以将其转换为均方根误差约为 2。这就是我的最终模型。
接下来,我可以问,嗯,这就是我能做到的最好结果吗?嗯,如果能有第二个模型进行比较就好了,这样我们就可以构建比我展示的这个模型复杂得多的模型。让我举一个更复杂的模型的例子,它看起来是这样的。之所以说它更复杂,是因为它有更多的参数,而这些参数又对应着一个更复杂的表达式。我不会深入讲解表达式的细节,我想重点关注拟合这个模型的机制。因为你会发现,如果你点击“运行”,它会失败。我想探究一下构建这些模型的一些原因和挑战。它提示我再次点击“运行”。再次点击“运行”,它仍然失败。这并没有真正帮助到我。
发生了什么?嗯,迭代次数已超过上限。上限是 60,这意味着什么?嗯,这意味着我们正在执行梯度下降算法。我就像坐在山坡上,正往下走,想到达山脚,但我只能走 60 步。嗯,我不太明白为什么会有这个限制。对于现代计算机来说,走超过 60 步没有任何问题。我们把限制设为 6000。我想重置这些参数。让我把它们恢复到初始值,然后再点击“运行”。它仍然失败,但我们得到了不同的错误信息。信息提示说这个收敛标准太严格了,我们需要放宽它。嗯,这个收敛标准的值非常小。我想 10 的负 15 次方不太合适。
我想要一个较小的数字,但我认为 10 的负 7 次方已经很小了。我们来试试,我再恢复到初始参数,然后点击“运行”。最终我收敛了,所以你可以看到,对于复杂的模型来说,有时需要更多的工作。现在我们来比较一下这两个模型。我的误差平方和是 128,或者说均方根误差大约是 2.5,而另一个模型的均方根误差是 2,误差平方和是 91。右边的模型比左边的模型表现更好。但至少我们可以比较不同的模型。不过,我们在进行比较时要谨慎,因为左边的模型仍然不是该模型的最佳表示。让我再给你展示一下我们还可以做些什么。
我要重新启动第二个模型的分析,但我不会点击“开始”按钮。我说过这是一个梯度下降算法,就像下山一样,我们讨论过可以走多少步,但这还取决于你的起点。如果你从不同的位置开始,最终会到达不同的位置。你的起点取决于这里的这个值。我想做的是探索不同的起始条件。JMP 有一个非常棒的功能,可以让我对这些起始条件进行网格搜索。我需要定义这个网格,但我不太擅长边打字边说话,所以我会打字。你可能会听到我说话,但我只是在自言自语,我需要设置好这个网格。
范围是 1 到 10,步长为 3。这是最大的,范围是 1 到 100,000,步长为 100。还有 0 到 10 的,我们取 10 个。好了,我絮絮叨叨说完了。现在我又跟你说话了。我已经设置好了。如果我点击“运行”,JMP 会评估所有这些初始条件。我有 27,000 个,我们来看看。嗯,不用全部都看,但你看,这里有 27,000 个。对于每一行,我都有这些参数的初始条件组合,它计算了误差平方和。现在我们把它们排序,按升序排列。你看,最佳误差平方和是 145。跟这个比起来还是不太好。但记住,我们还没有执行算法。这只是在寻找初始条件。
这些是我想要使用的初始参数。我复制一下,表格就完成了,关掉它。我可以把它们粘贴到这里。这是我的初始条件。现在我点击“运行”。现在我已经收敛了,我们来看看结果。可以看到,误差平方和是 1。记住,上次是 128,现在降到了 1。相比之前的 91,这简直太棒了。均方根误差是 0.2,而上次是 2,所以小了 10 倍。我在这里试图展示的是拟合这些模型的一些挑战,但其中蕴含着一个重要的教训:对于复杂的模型,如果不进行网格搜索来寻找初始条件,很可能会得到一个次优解。
我可以评估不同模型的各种参数,也可以比较不同的模型。我可能会使用均方根误差,但我也想对复杂模型进行惩罚,所以我至少可以在软件中计算AIC统计量,这是一个模型比较统计量。现在我要切换到软件,也就是JSL应用程序,它可以帮助我构建模型并进行模型比较。让我启动软件,界面是这样的。首先,我需要给一个新项目命名。就叫它“柏林”吧。我已经有一个“柏林”项目了。那就叫它“柏林2”吧,然后点击“确定”。你在这里看到的这些图标代表工作流程。首先,我们先来更全面地了解一下这个应用程序。
想想 JMP 的功能。JMP 使我们能够拟合模型、进行统计分析和可视化。当我们编写 JSL 脚本时,所有这些功能都相对容易实现,因为 JMP 已经为我们完成了所有繁重的工作。除此之外,此应用程序还旨在提供一种工作流程。我们需要执行一些步骤。这些步骤目前无法执行,它们已被锁定。这两个步骤是可用的,并且是推荐的步骤。我们提供了一个工作流程,即一系列要执行的任务,同时我们还会持久化信息。一旦设置完成,您就无需再次执行。如果您已经拟合了模型,则无需再次拟合。
你可以明天、下个月、明年再回来查看。这就是这个项目的理念:它能长期保存所有信息。我们先开始吧。第一步,也是推荐的一步,是选择产品。实际上,我在这里是在判断,我手里的产品是液体还是固体?用行话来说,液体产品是大分子产品,片剂是小分子产品。我的理解方式略有不同,我是从建模的角度来考虑的。液体对pH值很敏感,比如疫苗,它对pH值很敏感;而片剂则可能对湿度敏感。我打算吃一片片剂。我真希望我没说我要吃一片片剂。我真想笑出来。
接下来我要做的是分配一些数据,而我需要一些数据。这是一个数据表,它关联了四个杂质。回到应用程序。我可以分配这些数据,然后它会询问我是否要使用这些数据。是的,我同意。之后它会告诉我这些数据已经复制到项目中,我实际上不需要保持数据表打开状态。这些数据会作为项目的一部分持久化保存,但我们不需要实际查看该表。我们仍然可以进行列映射,所以我可以使用这四个杂质。这些是我的响应。这些是我要用来拟合模型的对象。对于每个杂质,我将拟合多个模型,然后我想为每个杂质选择最佳模型。
让我完成列映射。现在,对于每个响应,我都可以输入规格。我这里输入 1%,也就是 1%。我们想表达的是,如果杂质含量低于 1%,那么产品就是可行的。如果我们观察杂质含量随时间的变化,杂质含量超过 1% 的点就是产品的有效使用寿命。现在是最后一步,即提供模型选项。我们可以用八种不同的模型来拟合这些数据,这里有必要解释一下。想想杂质含量随时间的变化,它可能呈线性增长,也可能增长速度减慢,或者加快。我们有减速模型和加速模型的概念。
我们可以用不同的方式表达这些模型。我的鼠标好像丢了。好了,我们回到模型选项。这里有八种模型类型。这些模型包含温度和水分活度的表达式或项。水分活度可以用相对湿度或绝对湿度来表示。最终我们得到总共 16 种不同的模型变体。为了演示,我只使用八种模型变体,并且只关注涉及相对湿度的那些。我快完成了。还有一件事。数据本身是受压存储条件下的数据。我必须提供标准存储条件:25 摄氏度,65% 湿度。
项目设置任务到此完成。现在我可以开始拟合动力学模型了。应大家的要求,我可以进行自动拟合,只需点击确定即可。它会自动拟合所有模型并完成所有模型选择。但我希望你们真正理解模型选择的挑战。自动拟合并不能提供任何有价值的信息。我们至少要针对其中一个杂质单独进行模型选择。对于每种杂质,我们都要经历一个三步流程(或者我只以杂质一为例),找到合适的模型。第一步是拟合模型。程序会依次拟合每个模型。完成后,会显示一个仪表盘,确认这些模型已拟合完成。
最后一步是从这八个模型中选择一个。但最好先有一个中间步骤,创建一个候选名单。第二步就是创建这个候选名单。软件通过交通信号灯系统来简化这一步骤。绿色区域的模型是最有可能的模型,黄色区域的模型可能性较低,而红色区域的模型则可以舍弃。我们可以对这个名单进行合理性检查。例如,最佳模型是一个减速模型。如果杂质浓度随时间推移而降低,而不是加速,那么我们可以舍弃那些加速的模型。
软件并没有使用那种逻辑,也许它应该用,也许它应该具备那种智能。它现在只是冷冰冰的统计数据,计算AIC统计量。AIC是一种模型比较统计量,用于惩罚模型的复杂度。JMP在很多地方都用到了它,它贯穿整个JMP。虽然非线性平台本身并没有用到它,但我们可以在软件内部计算它。我接下来只关注绿色区域内的两个模型,我会把它们移除。我的候选名单就只有这两个模型了,现在我可以选择我最心仪的模型,或者至少软件会提供更多信息来帮助我做出最终选择。
我现在只关注这两个模型。现阶段统计数据本身对我帮助不大,因为我无法真正区分绿色区域内的两个模型。如果你看R平方值或均方根误差之类的指标,你会发现这些数值几乎完全相同。我能做的就是先把这些指标关掉,以免分散我们的注意力。专注于一些诊断图。我可以查看不同的诊断图,也可以查看我选择的两个模型的诊断图。我尝试在软件中加入一些提升用户体验的功能。我发现能够并排查看这两个模型非常有用。现在,我已经把第一个模型添加到列表中了。如果你愿意的话,最好的模型在左边,另一个模型在右边。
我们先明确一下我们关注的是什么:实际值与预测值。实际数据点与预测数据点的对比。如果我们有一个完美的模型,那么所有数据点都会落在45度红线上。这不仅能让你了解模型预测数据的效果,还能让我们比较两个模型。我希望你注意到的是,它们完全一样。它们如此相似。每次我看这些图的时候,我都怀疑是不是软件出了问题,没能给你展示不同的图表,但它们确实不一样。这个数据点在红线上,而这个数据点又在红线外。它们是不同的图表,但你可以看出它们有多么相似。这就是模型选择的挑战。这里很好地说明了这一点。
我们换个图看看。我们来看一下残差的正态概率图。我们希望残差能形成一条直线。在这个例子中,左边的图明显比右边的图好。左边的图也是统计结果最好的。这时我想,选择统计结果最好的模型,对吧?选择统计结果和诊断结果都最好的模型,你大概不会丢掉工作,所以没问题。但我还需要考虑如何使用这些模型。我可以看看预测分析器。我们知道预测结果是一样的,这一点我们已经看到了。这些分析器看起来会一样,或者至少在我们拥有的数据范围内看起来一样。但这并不是我们使用这些模型的方式。
首先,让我解释一下预测分析器,因为您可能希望看到温度和相对湿度。这些参数已移至下拉列表中。我们目前只关注时间演变。那么,我想做什么呢?我想将温度调回正常的存储条件,即 25 摄氏度和 65% 的湿度。然后,我想延长这些条件的时间。这又是一个提升用户体验的功能。我可以对这些图表进行外推。现在您可以看到显著的差异。当我们外推这些模型时,两个模型都描述了减速过程,但左侧的减速更为显著。
我们来看看这会带来什么后果。首先来看右边的模型。这是模型与这条线相交的点。这条线代表了该点的设定阈值。但这并不是我们关注的重点。当我们对寿命进行建模时,我们实际上关注的是置信上限以及它与模型的交点。这才是我们真正感兴趣的点。略高于 300 天。而左边的模型则更接近 600 天。这些数字非常重要,它们在软件中已经计算出来了,就在这里。但我希望你们理解它们的来源,以及它们在模型选择中的重要性。因为我们经常发现,这些是不同模型之间唯一有意义的统计差异。接下来我该怎么做呢?让我们思考一下我们为什么要进行建模。
我不是要告诉你这款产品的保质期。如果你想知道,把它放在温控箱里两年看看会发生什么。我们想做的是进行风险评估,看看这款产品是否有可能提前失效。这个模型告诉我,它有可能在300天后失效。当然,也可能不会。左边的模型可能是正确的,它的寿命可能达到600天。但我希望采取一种风险规避的方法。我可能会选择第二个模型,我可以解释原因。我可以说明,我采取的是基于风险的方法,而绿色区域的模型会给出一些解释。
如果你不满意,你可以说这是基于科学判断或其他原因。你也可以直接说我喜欢这个结果。但我们可以给出一些解释,或者我们可能想选择左边的模型,因为它统计数据最好。最终我们可以自己选择。但我认为这是基于风险评估的理性选择。我决定采用这个模型,这样就完成了对该杂质的模型选择。我们可以对所有杂质重复这个步骤,如果我给你看一个总结,你可以看到我们已经完成了杂质 1 的这三个步骤,我可以对其他三个杂质也做同样的事情。我可以手动操作,也可以使用自动拟合,现在这三个杂质都可以自动拟合了。
如果我这样做,我们可能要盯着屏幕看几分钟。不如这样,我来演示一下如何切换到另一个项目。我选择了另一个项目,数据是一样的,唯一的区别是所有四种杂质都进行了自动拟合。我只是想让你们看看输出结果。如果你们点击了自动拟合按钮,就会看到这样的输出结果。首先,杂质 1 的按钮是琥珀色的。软件里有很多类似交通信号灯的指示系统。其他杂质的按钮都是绿色的。琥珀色的按钮代表什么呢?它表示模型选择存在一些歧义,因为我可以选择这两个模型,而其他杂质的模型选择则没有歧义。
对于其他模型,我直接选择了统计值最小、统计值最佳的那个。而对于这个模型,我是根据持续完整性(CI)失效持续时间来选择的。这些持续时间都显示出来了。这是311,也就是我们之前看到的那个数字。如果我点击图标,或者只是把鼠标悬停在上面,它会告诉我估计值在311到569之间。如果我点击那里,它会显示交通灯系统,并告诉我它选择了第二个模型,因为它对这个数字的估计值最差。但情况并非总是如此。这里,这个模型的保质期非常短。我们看到一个绿色图标,绿色图标告诉我这里的估计值范围非常一致。绿色区域里有四个模型,但它们的预测值几乎相同。因此,我们可以根据统计数据选择最佳模型。
根据模型之间表现的一致性,我们会做出不同的决策。最终,我们完成了所有模型的拟合,并为每种杂质选择了合适的模型。这些杂质针对的是单一产品,之后我们可以进行一些分析。例如,我可以将这些模型汇总到仪表盘上,方便我们查看不同的模型。我们可以测试不同的储存条件。这是一种分析方法。另一种方法是将动力学模型应用于储存模型,并研究不同类型的容器对产品的影响。
让我总结一下我今天讲的内容。如果我没记错的话,我们首先来看一些来自加速稳定性研究的数据,这项研究旨在对产品施加压力。我提到的产品是一种药品,具体来说是药片。我可以模拟杂质的演变过程,构建这些动力学模型。我们已经了解了如何使用非线性平台生成这些模型,我也尝试向大家展示了使用该平台的一些挑战。然后,我讲解了如何使用 JSL 应用程序来自动化模型拟合并简化模型选择过程。你们还可以进一步思考包装产品的相关问题。今天的演讲就到这里。非常感谢。