Stata：经验动态建模

👇 连享会 · 推文导航 | www.lianxh.cn

• 🍎 Stata：Stata基础 | Stata绘图 | Stata程序 | Stata新命令
• 📘 论文：数据处理 | 结果输出 | 论文写作 | 数据分享
• 💹 计量：回归分析 | 交乘项-调节  | IV-GMM | 时间序列 | 面板数据  | 空间计量 | Probit-Logit | 分位数回归
• ⛳ 专题：SFA-DEA | 生存分析  | 爬虫 | 机器学习 | 文本分析
• 🔃 因果：DID | RDD | 因果推断 | 合成控制法 | PSM-Matching
• 🔨 工具：工具软件 | Markdown | Python-R-Stata
• 🎧 课程：最新专题 |  计量专题 | 关于连享会

🍓 课程推荐：2024 空间计量专题
主讲老师：范巧 (兰州大学)
课程时间：2024 年 10 月 2-4 日 (三天)
课程咨询：王老师 18903405450（微信）

作者：陈点点 (中国人民大学)
邮箱：chendiandian@ruc.edu.cn

编者按：本文主要摘译自下文，特此致谢！
Source：Jinjing L I, Zyphur M, Sugihara G. EDM: Stata module to implement empirical dynamic modeling[J]. 2019. -Link-

目录

• 1. 引言

• 2. 理论背景

• 2.1 单纯型投影

• 2.2 S - 映射

• 2.3 收敛交叉映射

• 3. Stata 实操

• 3.1 基于单纯形投影探索系统的维度

• 3.2 基于 S-map 的非线性检测

• 3.3 基于 CCM 的因果关系检测

• 4. 附录

• 5. 相关推文

温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：

1. 引言

如何研究不满足经典回归假设的复杂动态系统？在这些系统中，准自然实验或基于模型的回归往往无法适用，因果识别更是一项困难的任务。处理现实世界的复杂性需要能够描述和测试非线性动态系统因果关系的工具。

为了解决现有方法的局限性，研究人员开发了一种新方法，即经验动态建模 (empirical dynamic modeling, EDM)。这种来自自然科学的新方法允许：

• 描述动态系统，包括其复杂性、可预测性和非线性；
• 区分因果关系和单纯的相关性；
• 对非线性、稳定性或均衡所需假设最少。

Stata 中的 edm 命令可以实现这一方法。该命令为时间序列和面板数据提供了三种经验动态建模 (EDM) 方法：

• 单纯型投影，它能刻画系统的维度和确定性程度；
• S - 映射，用于量化系统的非线性程度；
• 收敛交叉映射 (CCM) ，它提供了一种非参数方法来建模因果效应。单纯型投影和 S - 映射通常以探索、诊断的方式用于描述产生观测时间序列或面板数据的系统，而 CCM 用于评估变量之间的因果关系。

接下来，我们将使用模拟数据来说明这些方法。

2. 理论背景

EDM 的逻辑基于这样一个事实，即产生数据的动态系统可以通过重构其基础系统来建模 (Takens 1981)。建议观看以下视频以直观了解 EDM 的思想：

• https://youtu.be/fevurdpiRYg
• https://youtu.be/QQwtrWBwxQg
• https://youtu.be/NrFdIz-D2yM

考虑一个具有 D 个时变变量的系统，例如一个沿着 GDP、失业和失业率变化的国民经济。随着时间的推移，D 个变量绘制了系统状态变化的轨迹图。下图是 D=3 情况的示意图 (呈现出 “蝴蝶” 吸引子)。随着国民经济发展，D 个变量的轨迹将在 D 维状态空间中形成 D 维流形 。流形 表示所有变量随时间变化的系统轨迹，因此在任何时间 t，系统状态都是 上反映系统变量 D 的单点。如果变量互为因果，则 代表一组典型的未知微分方程，这些方程生成一个 “吸引子”， 上的点将向该 “吸引子” 收敛。当然，吸引子可能是混沌的，而不是系统状态趋于收敛的固定点或点集 (均衡点)。其中，术语 “动态” 指的是以这种方式运行的系统。

如果吸引子流形 的基本方程已知，则基本上可以描述系统的复杂性和非线性，描述其确定性和随机性特征，并确定 D 变量之间的因果关系。然而，在实践中， 是未知的，并且所有 D 变量都无法测度。因此，必须使用时间序列或具有足够时间长度的面板数据中的观测变量 X 重构 。如果 X 是 的投影，Takens 定理表明，整个系统的确定性行为通常可以使用 X 的滞后来重构，从而形成 E 维阴影流形 (其中 $D<e<=2d+1$) x="" e="3">0$)。</e<=2d+1$)>

使用滞后来重构流形是一种状态空间重构的 “延迟嵌入” 或 “滞后坐标嵌入” 方法，其中 上的点由 E 长度的滞后向量 X 定义，重建质量通过样本外观测值和预测值之间的相关性 来评估。确定性系统的标志是预测。 反映了把基础系统重建为 的确定性流形的恢复程度。如果原始 D 维流形 在 E 空间中正确展开为 ，则预测能力 将最大化，因此跨越不同 E 值的  (即使用不同数量的滞后嵌入) 可用于推断基础系统 。这种用 重建 的一般方法适用大多情况，主要优点在于：

• 不需要 X 以外的额外变量来捕捉整个系统的动态；
• 对于确定性耦合的不可测变量，不需要进行控制；
• 未对线性、稳定性或均衡做出任何假设。

2.1 单纯型投影

为了将 重构为 ，使用 X 的 E 个滞后来构建一个长度为 E 的向量，该向量在 上形成一个单点 (即一个嵌入)，这是针对每个 进行的。在我们的方法中，默认情况下，长度为 E 的向量首先随机 50/50 分割，其中一半用于构建 的训练集。这样，训练集就是一个 E 长度的 X 滞后向量集，它的点落在重构的 E 维流形 上。

另一半数据形成 “预测” 或测试 / 验证集，其中包含落在 上某个位置的 E 长度 “目标” 向量。重构流形 中的信息可用于预测每个目标的未来。具体而言，对于预测集中的每个目标，训练集中 上的 最近邻 () 由欧几里德距离确定。这 k 个邻域在 上形成一个单纯形，用于将目标 包围在 E - 空间中。然后，包围目标的单纯形邻居被 “投影” 到未来 ，以计算预测目标未来值 的加权平均值。

每个邻居 i 的权重由目标到该邻居的欧几里德距离和距离衰减参数 确定。具体而言，权数 可以写为：

其中，

表示欧几里德距离， 是流形中最近的邻居 (即历史轨迹与目标最相似)。当 时，将忽略距离，并且所有邻居的权重相等。随着 的增加，附近邻居的权重增加，以更好表示局部状态 (即，与 X 更相似的历史轨迹)。默认情况下， 反映了更近邻居具有更大权重。

预测质量通过预测集的目标的未来实现值与投影单形的加权平均值的相关性 进行评估。 值越高，预测越准确。预测的平均绝对误差 (MAE) 是一种更侧重于观测数据和预测数据之间的绝对差距的度量，而不是像 这样的总体变化，也可以作为 的补充，具有相反的性质 (即 MAE 值越高，预测越不准确)。当变量被预先标准化时，MAE 范围将在 0 到 1 之间。

为了推断维度为 D 的基础系统，在不同的 E 值下评估 和 MAE，可以使用 –E 和 MAE–E 关系的函数形式来推断系统在研究时间内的确定性程度。与典型的回归方法不同，通过增加额外的滞后 (即，较大的 E) 来增加重构流形的维数，将增加影响预测的无关信息，降低预测准确度。因此， 在中等 E 时最大，并随着 E 的增加而下降。在具有自相关的高维或随机系统中，情况并非如此， 可能随 E 增加或随着 E 增加而接近渐近线，这就是为什么 E– 和 E–MAE 关系在诊断上是有用的。

理想情况下，一个系统可以用少于 10 个因素 (即，少于 10 个维度) 来描述，这样当 E<10 时预测最准确。在这种情况下，系统可被视为低维和确定性系统，其预测精度很高 (例如，>0.7 或 0.8)。

2.2 S - 映射

与单纯形投影一样，默认情况下，S - 映射使用 50/50 的数据分割为 E 长度嵌入向量的训练集和预测集。该训练集表示一个重构流形 ，流形上的 k 个最近邻用于预测预测集中每个目标向量的未来。对于 S - 映射，k 个相邻的训练集向量中的每个元素都可以被认为是类似于预测器的，考虑具有预测变量的 E 列的 k 行数据，使得预测器集合包括在 E 场合 。默认情况下包括常数项 c，这类似于具有 E+1 预测值和 k 观测值的局部回归，其中计算 E+1 系数以预测预测集中的每个目标。与单纯形投影不同的是，单纯形投影的邻域数为 k=E+1，而在 S - 映射中，通常选择 k 使其包括 上的整个点库 (即整个重构流形)。在数值上，点 t 处的预测值 y 计算公式为：

其中，系数向量 可以使用奇异值分解 (SVD) 计算，其形式为 B=AC，其中 B 是 k 维向量的加权未来值，所有相邻点和 A 是训练集中 k 个相邻点的权重矩阵 (包含原始时间序列中的过去和未来值，用于形成随机确定的训练集，用于将流形重建为 ) 以及常数项。数学上，，。S - 映射的权重 定义为：

当权重函数中的 时， 上的邻域没有不同的权重，因此在单变量情况下，S - 映射只是一个 E 阶自回归，训练数据与预测数据随机 50/50 分割  (即，映射是全局的，而不是局部的)。然而，随着权重 的增加，通过增加邻近邻居的权重进行预测，预测对系统的非线性行为变得更加敏感。换句话说，通过使用 上与预测集中目标更相似的历史轨迹中的更多信息，预测变得更依赖于状态。如果一个系统以依赖于状态的方式演化，那么来自附近邻居的更多信息将提高预测能力。

再次，使用 和 MAE，观察 – 和 MAE– 关系的函数形式，可以评估系统的性质。当 时，如果观察到状态依赖以较大的 和较小的 MAE 的形式存在，则可以使用 EDM 工具来建模系统的非线性动态行为。如果未观察到非线性，EDM 工具仍然可以使用 CCM 以非参数方式评估因果关系 (尽管在这种情况下，CCM 可能不如传统方法有效)。这里，S - 映射在诊断上可能是有用的，因为当 时，具有自相关的线性随机系统应显示最佳预测。如果没有其他原因，除了增加局部权重， 可能会增加对局部噪声的敏感性。

与单纯形投影一样，S 映射也可以以多变量方式进行。在多变量情况下，当 时，S - 映射更类似于自回归分布滞后 (ARDL) 或 DPD 模型。随着 的增加，它成为一种更局部的回归，预测越来越依赖于重构流形中的局部信息。在这里，标准化有助于确保模型中不同变量的权重相等，但应始终记住，S-map 不是传统回归模型，而是重建的吸引子流形 。

2.3 收敛交叉映射

收敛交叉映射基于这样一个事实，即如果 X 是 Y 的确定性驱动因素 (X→Y)，则 Y 的状态必须包含有助于恢复或 “交叉映射” X 状态的信息。该方法是单纯形投影的一个扩展，即使用一个变量重构吸引子流形 M，并使用它来预测不同的变量。如果变量共享同一个吸引子流形，则可以使用重构的流形进行预测。

具体来说，如果变量 X 和 Y 参与同一个动态系统流形 M，则重构的流形 和 可以相互映射。反过来，通过测试 X 和 Y 是否共享一个动态系统的信息，也可以测试变量在特定方向上相互影响的程度  (即 X→Y 和 / 或 Y→X) 。这是基于一个违反直觉的事实：如果 X→Y 是因果意义上的，那么关于 X 的历史信息包含在 Y 中，因此可以使用 通过单纯形投影预测 X，我们将其记为 。

这是违反直觉的，因为在典型的时间序列或面板数据方法中，原因是用来解释或预测结果，而不是相反。然而，在 CCM 中，结果变量 Y 交叉映射或 “xmaps” 原因变量 X，通过阴影流形 来预测 X  (即 ，可以望文生义地理解为潜在的 X→Y 影响) 。用结果变量预测原因变量，因为寻找原因需要从结果开始，并查看其动态结构 是否具有原因 X 的特征。注意，即使存在 X→Y 因果关系，如果 Y 不会引致 X，那么 在预测 Y 时的作用将很差，因为 将是 Y 以外变量的函数 (即， 效果不佳) 。

CCM 中的术语 “收敛” 描述了评估因果关系的标准。该术语反映了这样一个事实，即如果存在 X→Y 因果关系，则预测精度 (例如，X 和 之间的相关性 ) 将随着 上点的库大小 L 的增加而提高。在这种情况下，较大的库可以改进预测，因为它们使流形 M 更密集，最近邻变得更近。如果局部流形中存在因果信息，则可以改进预测。然而，如果 X–Y 关联仅仅是统计性的，那么增加 L 不会提高预测精度，因为更密集的流形不会提供额外的预测信息。总之，如果存在 X→Y 因果关系，那么 不仅会以 的预测精度的形式交叉映射 X，而且随着 上 L 点数量的增加，预测精度将提高。这一特性可以通过 L– (或 L–MAE) 关系进行检验。

3. Stata 实操

下面我们通过一个常见的动态系统来说明 EDM 命令的使用。两个变量 x 和 y 都依赖于各自的过去值，并通过方程的最后一项相互耦合。本例中，x 设置为仅由其自身的过去值确定，而 y 由 y 和 x 的过去值确定。

下图绘制了 和 分别设置为 0.2 和 0.3 时 x 和 y 随时间推移的值。删去 t 小于 300 的观测值，以使演化动态成熟。x 和 y 之间的相关系数为 0.15，在 0.05 水平上不显著，显示为两个不相关的变量。事实上，下图的曲线显示，在不同的时间窗口下，变量可能呈现正相关或负相关。这模拟了一种常见情况，即数据集中可能存在 “海市蜃楼相关性”：在不存在线性相关性的情况下，可能存在因果关系，而使用相关系数等典型线性指数可能会忽略这一点。

3.1 基于单纯形投影探索系统的维度

EDM 分析的第一步是确定系统的维度。维度可以理解为重建基础流形 所需的自变量数量。这是通过带有 explore 子命令的单纯形投影完成的，使用 e() 选项可以指定维数范围。

使用 explore 子命令，将数据随机分成两半，其中一半用作训练数据集以构造阴影流形 ，另一半是用于评估 上投影的样本外预测能力的测试数据集。通常根据预测值和观测值之间的最高 或最低 MAE 来选择最佳 E  (以保持模型的简约性) 。这两个结果通常一致，如果不一致，则使用表示最低嵌入维数的 E。此外，在小样本的情况下，MAE 是首选，因为它对异常值不太敏感。

在本例中，我们使用单纯形投影探索 2 到 10 之间的所有维度，以确定最佳 E。对于大型数据集，可以选择 E=20 或更高。由于训练数据集是随机选择的，我们将该方法复制了 50 次，通过选项 rep(50) 来估计 50 次运行的平均性能。命令和输出如下所示。ci() 选项可用于报告从多次复制中得出的 平均值的估计置信区间。edm 命令也可与 jackknife pre fix 结合使用，以实现附加控制。通常选择具有最高 或最低 MAE 的 E 值，因为其重建的流形与观测数据最匹配。

. edm explore y, e(2/10) rep(50)
Replication progress (50 in total)
.................................................. 50

Empirical Dynamic Modelling
Univariate mapping with y and its lag values
----------------------------------------------------------------------
                      --------- rho --------- --------- MAE ---------
 Actual E    theta         Mean    Std. Dev.         Mean    Std. Dev.
----------------------------------------------------------------------
        2        1       .97999     .0059908      .031708      .004511
        3        1       .97083     .0078036      .039437     .0044651
        4        1       .95058      .013928      .049484     .0060995
        5        1       .92494      .023275      .059684     .0078145
        6        1       .89361      .032127      .069651     .0094557
        7        1        .8451      .046727      .082331      .010843
        8        1       .79302      .052859      .096854      .010539
        9        1       .75721      .058412       .10613      .010663
       10        1       .72319      .063342       .11468      .010283
----------------------------------------------------------------------
Note: Results from 50 runs
Note: Number of neighbours (k) is set to E+1
Note: Random 50/50 split for training and validation data

结果表明，随着 E 的增加， 下降，MAE 增加。这表明最优 E 为 2，正好是我们用来创建动态系统的自变量数量。还可以使用返回的 e(explore_result) 矩阵的内容绘制结果，如下图所示。

3.2 基于 S-map 的非线性检测

如果流形 M 上的轨迹相对于系统当前状态不变，则存在线性；而如果系统演化依赖于状态，则存在非线性。通过从单纯形投影中选择最佳 E 后，接下来我们可以估计一个回归模型，该模型使用距离衰减参数 改变附近观测值的权重 (根据系统状态而不是时间)。然而，我们再次强调，不应将其解释为典型的回归模型，而应将其解释为重建的流形。

当选项 theta(0) 中 时， 上的邻域没有不同的权重，因此 S-map 简化为一种自回归模型，训练数据与测试数据随机 50/50 分割。然而，随着 的增加，预测对系统的非线性行为变得更加敏感，因为它更多地利用附近的观测数据进行预测。换句话说，预测变得更加依赖于状态。再次使用 和 MAE 并形成 和 MAE– 图，可以评估系统的非线性。当 时，如果以较大 和较小 MAE 的形式观察到非线性或状态依赖性，则需要 EDM 工具来建模系统行为。如果系统是线性的， 不会随着 超过零而增加。在这种情况下，具有线性假设的模型也可能是合适的。

explore 和 xmap 子命令都支持 S 映射，而不是局部预测的单纯形投影。下面给出了一个例子来分析前一个例子中变量 y 的非线性。在这种情况下，我们探索 0 和 5 之间可能的 值，增量为 0.01。此外，我们还通过在 k() 选项中指定一个负数来包含本地预测的所有观测值。这允许在低 E 或 的低数据密度区域获得更稳定的结果。

. edm explore y, e(2) algorithm(smap) theta(0(0.01)5) k(-1) 

Empirical Dynamic Modelling
Univariate mapping with y and its lag values
--------------------------------------------------------------------
          Actual E           theta             rho             MAE
--------------------------------------------------------------------
                 2               0          .67288          .13505 
                 2             .01          .67831          .13412 
                 2             .02          .68365           .1332 
                 2             .03           .6889          .13228 
(output omitted)
                 2            4.94          .99378         .019535 
                 2            4.95          .99379         .019527 
                 2            4.96          .99379          .01952 
                 2            4.97           .9938         .019513 
                 2            4.98          .99381         .019505 
                 2            4.99          .99382         .019498 
                 2               5          .99382          .01949 
--------------------------------------------------------------------
Note: Number of neighbours (k) is adjusted to 74
Note: Random 50/50 split for training and validation data

上面命令的结果也可以使用返回矩阵 e(explore_result) 以图形方式绘制，如下图所示。 随 的增加而增加，这与系统高度非线性的特性是一致的。

3.3 基于 CCM 的因果关系检测

在前面的步骤中，分析表明本例的最佳 E 为 2，与数据生成过程一致。现在，我们使用 xmap 子命令，使用之前观察到的 e=2 值，推导变量 x 和 y 之间的双变量 (整体) 因果效应。虽然在大多数情况下，我们建议将变量标准化，使其处于相同的尺度，但在这个假设的示例中，标准化对结果影响不大。

CCM 使用一个变量的重构流形 (潜在结果) 预测另一个变量 (潜在原因) 的值来评估变量之间的偶然联系。本质上，这是一种双变量单纯形投影，使用一个变量形成训练集，同时对另一个变量进行预测。使用预测值与观测值，计算 和 MAE。双变量交叉映射捕获了总体偶然效应，包括直接效应 (X→Y) 和可能的间接效应 (X→Z→Y)。

值得注意的是，单纯形投影结果可能表明不同的变量应使用不同的 E 值。如果是这种情况，选择用于重建流形的 E 可与 direction(oneway) 选项一起使用，以对每个变量运行单独的 CCM 程序。在这种情况下，应该采用原因变量的 E，因为原因变量是由结果预测的。例如，如果对 y 而言 E=2，对 x 而言 E=4，则可以运行以下 CCM：edm xmap x y，e (2) direction (oneway)，用于 或 的情况；而 edm xmap y x，e (4) direction (oneway) 用于 或 情况。由第一个变量重构的流形用于预测第二个变量，单纯形投影的 E 适用于所列的第二个变量。

. edm xmap x y, e(2)

Empirical Dynamic Modelling
Convergent Cross-mapping result for variables x and y
---------------------------------------------------------------
      Mapping    Library size             rho             MAE
---------------------------------------------------------------
   y ~ y|M(x)             150          .23019          .19673 
   x ~ x|M(y)             150          .69682          .13714 
---------------------------------------------------------------
Note: The embedding dimension E is 2

默认情况下，在不为 library() 指定值的情况下，整个数据集将用作训练库。换句话说，来自一个变量的所有可能信息都用于预测另一个变量。这通常是 CCM 的第一步，因为事先不容易知道最优 L 是什么。此外，如果最大 L 的 CCM 没有显示出有意义的可预测性，那么收敛的问题就没有意义了。然后，可以在评估整个结果时手动设定最大 L。具体而言，如果存在因果关系 (即 ) ，我们期望在最大 L 时进行有意义的预测，并且随着库中使用更多数据点，预测应得到改善。为了评估这一点，library() 选项指定了库大小的范围，replicate() 选项允许重复估计。当指定的 L 小于最大可用库长时，将为重建流形选择一个随机子样本，并且可以选择 replicate() 选项。

下面的示例估计了库长为 5 到 150 时的 ，并重复该过程 10 次 (在每个库长下随机抽取 10 次以形成流形) 。多次重复是为了确保有足够的点可用于估计，以清楚地识别所有趋势，以及映射的两个方向 (x–y 和 y–x；见下图) 之间趋势的差异。当然，对于最大库长，因为使用了所有观察值，此复制过程不会产生不同的结果。最大库长在计算上要求最高，所以对于大型数据集，谨慎的做法是仅在最大库长下单独运行一次 xmap 过程，然后在较小的库长下使用复制选项来评估不确定性。

. edm xmap x y, e(2) rep(10) library(5/150)
Replication progress (20 in total)
....................

Empirical Dynamic Modelling
Convergent Cross-mapping result for variables x and y
---------------------------------------------------------------
      Mapping        Lib size        Mean rho       Std. Dev.
---------------------------------------------------------------
   y ~ y|M(x)               5       -.0090982          .12325
   y ~ y|M(x)               6        -.004514          .11706
   y ~ y|M(x)               7         .018604          .10854
   y ~ y|M(x)               8          .01794          .10815
   y ~ y|M(x)               9         .033454         .095494
   y ~ y|M(x)              10         .037815         .097865
   y ~ y|M(x)              11         .030339          .10197
   y ~ y|M(x)              12         .040244          .10324
   y ~ y|M(x)              13         .049327         .098762
(output omitted)
   x ~ x|M(y)             142          .68351         .010373
   x ~ x|M(y)             143          .68683        .0093171
   x ~ x|M(y)             144          .68719        .0092867
   x ~ x|M(y)             145           .6898         .011021
   x ~ x|M(y)             146          .68971        .0098815
   x ~ x|M(y)             147           .6911        .0073416
   x ~ x|M(y)             148          .69069        .0061557
   x ~ x|M(y)             149          .69297        .0054071
   x ~ x|M(y)             150          .69682               0
---------------------------------------------------------------
Note: Results from 10 replications
Note: The embedding dimension E is 2

详细的结果存储在两个返回矩阵 e(xmap_1) 和 e(xmap_2) 中，可以使用 Stata 的 svmat 命令将其转换为变量，以便于操作。e(xmap_1) 将变量 x 视为变量 y 的结果 (即 ) ，e(xmap_2) 将变量 y 视为变量 x (即 ) 的结果。下面的示例存储结果并绘制 –L 收敛图，如下图所示。每个绘制的点代表一个点估计值，每个库长有十个估计值。局部多项式平滑线显示了 值随重建流形变得更密集 (即 L 增加) 的一般趋势。

–L 诊断图表明，从 y (即 ) 构建的流形预测的 x 始终优于 ，表明数据序列 y 包含更多关于 x 的信息。这表明 x 是 y 的 CCM 原因 (即 X→Y) ，与模型的实际情况相一致：y 由 x 引起，x 完全由其滞后值而不是 y 决定。

还应注意，当两种预测相似且处于预测性能都较好时，结果不排除双向因果关系的可能性。在这种情况下，预测性能的相对差异可用于确定更主要的因果方向，而不排除双向因果关系的推断。当从两个方向的预测性能都较差时，则变量可能是因果独立的。当两个变量的可预测性都很高，但没有观察到收敛，那么有可能存在一个外生的第三个变量可在驱动这两个变量。

4. 附录

/* edm paper replication do file
Jinjing Li, University of Canberra
Michael Zyphur, University of Melbourne
*/

clear all
set more off
set scheme sj

/* Example with a Synthetic Dataset */
/* Create a dynamic system */

set obs 500
gen t = _n
tsset t
gen x = 0.2 if _n==1
gen y = 0.3 if _n==1
gen z = 0.1 if _n==1
local r_x 3.79
local r_y 3.79
local r_z 3.77
local beta_xy = 0.0
local beta_yx=0.2
local beta_zy = 0.5
local tau = 1
drawnorm u1 u2
qui {
    forvalues i=2/`=_N' {
        replace x=l.x *(`r_x' *(1-l.x)-`beta_xy'*l.y) in `i'
        replace y=l.y *(`r_y' *(1-l.y)-`beta_yx'*l`tau'.x) in `i'
        replace z=l.z *(`r_z' *(1-l.z)-`beta_zy'*l`tau'.y) in `i'
    }
}
keep in 300/450
set seed 12345678

/* Typical Output */
edm explore x
edm xmap x y

/* Plot the system */
twoway (connected x t)(connected y t) in 1/40, ytitle(value) xtitle(t) ///
       legend(position(7) ring(0)) 
graph export plot.pdf, replace

/* rho-E plot */
edm explore y, e(2/10) rep(50)
mat r= e(explore_result)
svmat r, names(col)
twoway (scatter c3 c1)(lpoly c3 c1),xtitle("E") ytitle("{it:{&rho}}") ///
       legend(order(1 "{it:{&rho}}" 2 "local polynomial smoothing")   ///
       col(1) position(8) ring(0))
graph export rho-e.pdf, replace
drop c*

/* rho-theta plot */
edm explore y, e(2) algorithm(smap) theta(0(0.01)5) k(-1) 
mat r= e(explore_result)
svmat r, names(col)
twoway (line c3 c2) , legend(order(1 "{it:{&rho}}") position(5) ring(0)) ///
       xtitle("{it:{&theta}}") ytitle("{it:{&rho}}") ///
       title("{it:{&rho}}-{it:{&theta}} of variable y")
graph export rho-theta.pdf, replace
drop c*

/* ccm */
edm xmap x y, e(2)

/* ccm-convergent plot */
edm xmap x y, e(2) rep(10) library(5/150)
mat c1= e(xmap_1)
mat c2= e(xmap_2)
svmat c1, names(xy)
svmat c2, names(yx)
label variable xy3 "y|M(x)"
label variable yx3 "x|M(y)"
twoway (scatter xy3 xy2, mfcolor(%30) mlcolor(%30)) ///
    (scatter yx3 yx2, mfcolor(%30) mlcolor(%30))    ///
    (lpoly xy3 xy2)(lpoly yx3 yx2), xtitle(L) ytitle("{it:{&rho}}") legend(col(2))
graph export rho-L.pdf, replace

5. 相关推文

Note：产生如下推文列表的 Stata 命令为：
  lianxh 动态 系统 edm, m
安装最新版 lianxh 命令：
  ssc install lianxh, replace

• 专题：计量专题
• 知乎热议：如何系统地自学计量经济学？
• 专题：Stata入门
• 珠联璧合 I：Jupyter Notebook 和 Stata 关联 (windows系统)
• 专题：Stata绘图
• Stata绘图：制作教学演示动态图-GIF
• 专题：结果输出
• Stata 可重复性报告系列A：动态文档命令 (dyn*)
• 专题：面板数据
• 面板数据模型-xtdcce2：动态共同相关和截面相关
• Stata：动态面板数据模型OLS估计的偏差
• Stata实操陷阱：动态面板数据模型
• Stata: 动态面板门槛模型
• 专题：时间序列
• DCC-MGARCH：动态条件相关系数模型 (R+Stata)
• 专题：空间计量
• 刘迪：Stata空间溢出效应的动态图形-空间计量
• 空间计量溢出效应的动态GIF演示
• 专题：Probit-Logit
• 动态 Probit 模型及 Stata 实现