计算系统发育数据的TIGER率

这是一个纯Python实现的TIGER率（Cummins和McInerney 2011），可以用于计算所有类型的系统发育数据。该包的创建灵感来源于Syrjänen等人（2021）的一项研究，其中作者提出了一个用于语言目的的TIGER率计算的新实现。由于Syrjänen等人的实现不是一个包，而是一个命令行工具，其从Python脚本中的应用并不直接，因此我决定根据两篇论文中的描述重新实现率计算，以一个小型Python包的形式，不包含其他依赖项。

工作流程

为了计算您的系统发育特征的TIGER率，您的特征应以嵌套Python字典的形式存储，其中每个特征由一个键表示，该键链接到一个表示特征状态的字典。特征状态本身由一个以分类单元（语言变体、物种）为键，特征状态作为列表形式值表示的字典表示。这允许每个分类单元具有多态特征，并允许缺失数据（作为空列表传递）。在这种形式下，特征状态可以方便地以JSON格式存储，例如以下来自Cummins和McInerney的示例：A = CTTAA，B = AGGGG

{
  "A": {
    "1": ["C"],
    "2": ["T"],
    "3": ["T"],
    "4": ["A"],
    "5": ["A"]
  },
  "B": {
    "1": ["A"],
    "2": ["G"],
    "3": ["G"],
    "4": ["G"],
    "5": ["G"]
  }
}

将字符以这种形式存储后，您可以将其转换为集合划分

>>> from pylotiger import get_set_partitions
>>> characters = {"A": {"1": ["C"], "2": ["T"], "3": ["T"], "4": ["A"], "5": ["A"]}, "B": {"1": ["A"], "2": ["G"], "3": ["G"], "4": ["G"], "5": ["G"]}}
>>> taxa = ["1", "2", "3", "4", "5"]
>>> set_partitions = get_set_partitions(characters, taxa)

输出是一个字典，字符作为键，相应的集合划分作为值。集合划分由包含分类单元（在此以数字1到5表示）的 frozensets 组成。

>>> set_partitions
{'A': {frozenset({'1'}), frozenset({'4', '5'}), frozenset({'2', '3'})},
 'B': {frozenset({'1'}), frozenset({'2', '3', '4', '5'})}}

由于我们的实现允许多态字符，不同的集合划分可能包含相同的分类单元。如果某个字符的状态信息缺失，相应的分类单元将不会出现在任何集合划分中。

计算了集合划分后，我们可以将它们相互比较。Cummins 和 McInerney 报道，比较A与B得到得分为0.5，而比较B与A得到得分为1。

>>> from pylotiger import get_partition_agreement_score
>>> get_partition_agreement_score(set_partitions["A"], set_partitions["B"])
0.5
>>> get_partition_agreement_score(set_partitions["B"], set_partitions["A"])
1

为了计算速率，我们现在使用 Syrjänen 等人（2021年）研究的示例数据。

{
  "X": {
    "a": ["A"],
    "b": ["B"],
    "c": ["A"],
    "d": ["A"],
    "e": ["C"]
  },
  "Y": {
    "a": ["D"],
    "b": ["D"],
    "c": ["E"],
    "d": ["E"],
    "e": ["F"]
  },
  "Z": {
    "a": ["G"],
    "b": ["H"],
    "c": ["H"],
    "d": ["H"],
    "e": ["H"]
  }
}

在这里，作者报告了特征X和Y的速率为0.5835，特征Z的速率为0.667。

>>> from pylotiger import get_rates
>>> characters = {"X": {"a": ["A"], "b": ["B"], "c": ["A"], "d": ["A"], "e": ["C"] }, "Y": {"a": ["D"], "b": ["D"], "c": ["E"], "d": ["E"], "e": ["F"] }, "Z": {"a": ["G"], "b": ["H"], "c": ["H"], "d": ["H"], "e": ["H"]}}
>>> taxa = ["a", "b", "c", "d", "e"]
>>> set_partitions = get_set_partitions(characters, taxa)
>>> rates = get_rates(set_partitions)
>>> for char in "XYZ":
...     print(char, "{0:.4f}".format(rates[char]))
X 0.5833
Y 0.5833
Z 0.6667

请注意，在此版本中，速率计算没有针对重复字符值进行优化。因此，当使用它来研究具有许多字符的数据集时，当前实现可能需要一些时间。

计算修正的TIGER速率

TIGER速率旨在测试给定数据集的树状特征。在Cummins和McInerney提出的形式中，不变字符（在整个分类单元中具有相同字符状态）并未明确排除在数据之外。这可能是问题，因为不变字符不提供关于子分组的具体信息。然而，在原始划分一致性得分中，它们在与任何其他集合划分比较时总会得到得分为1，正如我们在前一个例子中将我们的字符“X”设置为反映相同的字符状态并与其他字符比较时可以容易地看到的那样。

>>> characters = {"X": {"a": ["A"], "b": ["A"], "c": ["A"], "d": ["A"], "e": ["A"] }, "Y": {"a": ["D"], "b": ["D"], "c": ["E"], "d": ["E"], "e": ["F"] }, "Z": {"a": ["G"], "b": ["H"], "c": ["H"], "d": ["H"], "e": ["H"]}}
>>> taxa = ["a", "b", "c", "d", "e"]
>>> set_partitions = get_set_partitions(characters, taxa)
>>> for char in "YZ":
...     print("X vs.", char, "{0:.4f}".format(get_partition_agreement_score(set_partitions["X"], set_partitions[char])))
X vs. Y 1.0000
X vs. Z 1.0000

相反，如果一个字符对于所有分类单元都有不同的状态，所有其他的数据集中 其他 集合划分在与之比较时都将得到得分为1，正如我们可以通过修改前一个例子中的字符“X”使其对所有分类单元具有不同的状态时容易地看到的那样。

>>> characters = {"X": {"a": ["A"], "b": ["B"], "c": ["C"], "d": ["D"], "e": ["E"] }, "Y": {"a": ["D"], "b": ["D"], "c": ["E"], "d": ["E"], "e": ["F"] }, "Z": {"a": ["G"], "b": ["H"], "c": ["H"], "d": ["H"], "e": ["H"]}}
>>> taxa = ["a", "b", "c", "d", "e"]
>>> set_partitions = get_set_partitions(characters, taxa)
>>> for char in "YZ":
...     print(char, "vs. X", "{0:.4f}".format(get_partition_agreement_score(set_partitions[char], set_partitions["X"])))
Y vs. X 1.0000
Z vs. X 1.0000

有人可能会认为这种行为与原始的速率计算一致，但当我处理具体的语言数据时，我发现包含对于子分组不提供信息的字符有时会扭曲得分，给人一种数据集高度树状的印象，仅仅因为它包含许多具有单例字符状态的字符。因此，此TIGER速率实现提供了划分一致性得分的修正版本，它在与原始划分一致性得分在三个方面有所不同

1. 单例（在数据中只出现一次的字符状态）被排除在兼容性计算之外，
2. 不变字符状态（对所有字符相同）被排除在计算之外，
3. 而不是将匹配项（其中划分A包含划分B）除以划分B中的集合数，我们将其除以交集数，这允许捕捉字符状态与两个划分兼容的情况。

关于1和2的处理，我们有两种可能性：如果我们遇到一个只包含单例状态或只具有不变状态的字符，我们可以将得分设置为0，或者我们可以从计算中忽略它。

要计算修正的TIGER得分，只需要导入 corrected_pas() 函数，并将其作为 partition_func 参数传递给 get_rates() 函数。

>>> from pylotiger import corrected_pas
>>> characters = {"X": {"a": ["A"], "b": ["B"], "c": ["C"], "d": ["D"], "e": ["E"] }, "Y": {"a": ["D"], "b": ["D"], "c": ["E"], "d": ["E"], "e": ["E"] }, "Z": {"a": ["G"], "b": ["G"], "c": ["G"], "d": ["H"], "e": ["H"]}}
>>> taxa = ["a", "b", "c", "d", "e"]
>>> set_partitions = get_set_partitions(characters, taxa)
>>> for x in "XYZ":
...     for y in "XYZ":
...         pas1 = get_partition_agreement_score(set_partitions[x], set_partitions[y])
...         pas2 = corrected_pas(set_partitions[x], set_partitions[y])
...         print("{0} | {1} | {2:.2f} | {3}".format(x, y, pas1, "{0:.2f}".format(pas2) if pas2 is not None else "None"))
X | X | 1.00 | None
X | Y | 0.00 | None
X | Z | 0.00 | None
Y | X | 1.00 | None
Y | Y | 1.00 | 1.00
Y | Z | 0.50 | 0.33
Z | X | 1.00 | None
Z | Y | 0.50 | 0.33
Z | Z | 1.00 | 1.00

从本例中可以看出，我们排除了所有与字符X的比较，因为该函数返回None。在计算TIGER比率时，这些字符将不会包含在比率计算中。此外，计算字符Y（状态：DDEEE）和Z（状态：GGGHH）的得分也显示出关键性的差异，因为我们在Z的字符状态和Y的字符状态之间识别出三种可能的链接（G=abc → D=ab，G=abc → E=cde，和H=de → E=cde），其中只有一种是兼容的（G → D）。因此，得分被设置为0.33（1/3），而不是原始TIGER分区协议比率计算中的0.5。