Julia 0.7 中一流的统计缺失值支持

Julia 的 0.7 版本将很快引入对统计缺失值的一流支持。作为统计分析和数据管理必不可少的功能，此功能在诸如SQL（使用NULL值）、R（NA）、SAS（.、' '等）或Stata（.等）等专门语言中很常见。然而，它在通用语言中却很少见，在通用语言中，Nullable或Option类型通常不允许隐式传播空值（它们需要提升），并且不提供具有缺失值的数组的有效表示^[1]。

从 Julia 0.7 开始，缺失值使用新的missing对象表示。它是经过多年来持续的设计讨论、实验和语言改进的结果，是DataArrays包中实现的NA值的继承者，该包过去是 Julia 中表示缺失数据的标准方式。missing实际上与其前身NA非常相似，但它受益于 Julia 编译器和语言的许多改进，使其速度更快，从而可以放弃DataArray类型并改用标准Array类型^[2]。借鉴现有语言的经验，missing的设计紧随 SQL 的NULL和 R 的NA的设计，后者可以被认为是在缺失值支持方面最一致的实现。顺便说一下，这使得从 Julia 代码生成 SQL 请求或使 R 和 Julia 互操作变得容易。

此框架由版本 0.11的DataFrames包使用，该包即使在性能改进仅在 Julia 0.7 上可用时，也已经通过Missings包在 Julia 0.6 上运行。

统计缺失值的新实现遵循我们认为对于现代科学计算环境至关重要的三个原则

缺失值默认情况下是安全的：当传递给大多数函数时，它们要么传播要么抛出错误。
missing对象可以与任何类型结合使用，无论是定义在 Base 中、包中还是用户代码中。
使用缺失值的标准 Julia 代码效率很高，无需特殊技巧。

这篇文章首先介绍了新的missing对象的行为，然后详细介绍了它的实现，特别是展示了它如何在保持完全通用性的同时提供极快的性能。最后，讨论了当前的限制和未来的改进。

“missing”对象：安全且通用的缺失值

Julia 的优势之一是用户定义的类型与内置类型一样强大且快速。为了充分利用这一点，缺失值不仅必须支持Int、Float64和String等标准类型，还必须支持任何自定义类型。出于这个原因，Julia 无法像 R 和 Pandas 一样使用所谓的哨兵方法来表示缺失，即在类型的域内保留特殊值。例如，R 使用可表示的最小 32 位整数（-2,147,483,648）来表示整数和布尔向量中的缺失值，并使用特定的NaN有效负载（1954，据说指的是 Ross Ihaka 的出生年份）来表示浮点向量中的缺失值。Pandas 仅支持浮点向量中的缺失值，并将它们与NaN值混淆。

为了提供可以与任何类型组合的缺失值的一致表示，Julia 0.7 将使用missing，一个没有字段的对象，它是Missing单例类型的唯一实例。这是一种普通的 Julia 类型，为此实现了一系列有用的方法。可以是T类型或缺失的值可以简单地声明为Union{Missing,T}。例如，保存整数或缺失值的向量类型为Array{Union{Missing,Int},1}

julia> [1, missing]
2-element Array{Union{Missing, Int64},1}:
 1
  missing

这种方法的一个有趣的特性是，一旦缺失值被替换或跳过（见下文），Array{Union{Missing,T}}的行为就像普通的Array{T}一样。

如上例所示，定义了提升规则，以便连接类型为T的值和缺失值会生成一个元素类型为Union{Missing,T}而不是Any的数组。

julia> promote_type(Int, Missing)
Union{Missing, Int64}

这些提升规则对于性能至关重要，我们将在下面看到。

除了通用和高效之外，新missing框架的主要设计目标是确保安全性，即缺失值永远不会被静默忽略或替换为非缺失值。缺失值是统计工作中的一个棘手问题，并且是错误或无效结果的常见来源。忽略缺失值相当于执行数据插补，这永远不应该在没有明确请求的情况下静默发生。不幸的是，某些主要的统计语言就是这样：例如，在 SAS 和 Stata 中，即使x缺失，x < 100也会静默返回true或false^[3]。众所周知，这种行为会导致已发表的科学著作中的结果不正确^[4]。哨兵方法也存在一些极端情况下的错误：例如，在 R 中，NA + NaN返回NA，但NaN + NA由于浮点计算规则返回NaN。

因此，将missing传递给函数将始终返回missing或抛出错误（除了下面介绍的一些特殊函数）。为了方便起见，标准运算符和数学函数系统地传播缺失值。

julia> 1 + missing
missing

julia> missing^2
missing

julia> cos(missing)
missing

julia> round(missing)
missing

julia> "a" * missing
missing

归约操作继承了基本运算符的传播行为。

julia> sum([1, missing, 2])
missing

julia> mean([1, missing, 2])
missing

另一方面，使用缺失值索引到Vector中是一个错误。缺失值不会被静默跳过，这相当于假设它们为false。

julia> x = 1:3
1:3

julia> x[[true, missing, false]]
ERROR: ArgumentError: unable to check bounds for indices of type Missing

julia> x[[1, missing]]
ERROR: ArgumentError: unable to check bounds for indices of type Missing

提供了便利函数来明确地去除缺失值。首先，skipmissing函数返回传递的集合中非缺失值的迭代器。在计算归约时忽略缺失值时，它特别有用。调用collect以获得包含所有非缺失值的向量。

julia> sum(skipmissing([1, missing, 2]))
3

julia> mean(skipmissing([1, missing, 2]))
1.5

julia> collect(skipmissing([1, missing, 2]))
2-element Array{Int64,1}:
1
2

julia> x[collect(skipmissing([1, missing, 2]))]
2-element Array{Int64,1}:
1
2

其次，coalesce函数返回第一个非缺失参数（如 SQL 中），作为特例，它允许用特定值替换缺失值。结合“点”广播语法，它允许替换数组中的所有缺失值。

julia> coalesce(missing, 0)
0

julia> coalesce(missing, missing, 0)
0

julia> coalesce.([1, missing, 2], 0)
3-element Array{Int64,1}:
1
0
2

julia> coalesce.([1, missing, 2], [2, 3, missing])
3-element Array{Int64,1}:
1
3
2

一组受限的函数和运算符遵循与上面描述的不同语义。它们可以分为四类。

ismissing如果输入为missing，则返回true，否则返回false。
===、isequal和isless始终返回布尔值。===和isequal在将missing与missing进行比较时返回true，否则返回false。isless也属于此类，并在将任何非缺失值与missing进行比较时返回true，否则返回false：缺失值最后排序。
&、|和⊻/xor实现了三值逻辑，根据结果是否即使在不知道missing背后的值是什么的情况下也能完全确定，返回布尔值或missing。
==、<、>、<=和>=如果其中一个操作数为missing，则返回missing，就像任何其他运算符一样。当在包含缺失值的集合上调用时，这些运算符会递归应用并遵循三值逻辑：如果结果取决于missing元素将取什么值，则返回missing。这也适用于all、any和in。

短路运算符&&和||，就像if条件一样，如果需要评估缺失值，则会抛出错误：在这种情况下，无法确定是否应运行代码。

当然，基于上述列出内容定义的函数继承了它们的行为。例如，findall(isequal(1), [1, missing, 2])返回[1]，但findall(==(1), [1, missing, 2])在遇到缺失值时会抛出错误。

有关这些规则的更多详细信息和说明，请参阅手册。如上所述，它们通常与 SQL 的NULL和 R 的NA实现的规则一致。

从“Nullable”到“missing”和“nothing”

虽然它类似于以前的NA值，但新的missing对象也替换了 Julia 0.4 中引入的Nullable类型，该类型被证明不是表示缺失值的最佳选择^[5]。Nullable存在几个问题。

它用于表示两种截然不同的缺失类型，我们有时分别称之为“软件工程师的空值”和“数据科学家的空值”。前者指的是Nullable在大多数语言中通常用于表示什么，即null表示值不存在，并且Nullable类型的一个优点是强制开发人员在没有值时显式处理这种情况。后者指的是统计缺失值，这些值通常在专门语言中静默传播，这对于方便性至关重要。越来越清楚的是，这两种用法存在冲突（即使专门的语法可以帮助缓解这种情况）。
为了确保类型稳定性，无论值是否存在，都必须指定Nullable{T}的T类型参数。在没有任何值的情况下找出合适的类型（称为“反事实返回类型”）在许多情况下已被证明是有问题的，在这些情况下，代码不得不严重依赖对类型推断的显式调用，这些调用最好由编译器直接处理。
Array{Nullable{T}}对象使用次优的内存布局，其中T值和关联的Bool指示符并排存储，这由于对齐约束而浪费空间，并且对于处理而言效率不高。因此，必须使用诸如NullableArray之类的专门数组类型（类似于DataArray）。

由于所有这些原因，Nullable在 Julia 0.7 中不再存在。根据用例提供了多个替换项。

如上所述，可能包含统计缺失值的数据应表示为Union{Missing,T}，即T类型的值或missing对象。
可能存在T类型的值或不存在值的情况应使用Union{Nothing,T}（在 Julia 0.6 中等效于Union{Void,T}）。作为特例，如果nothing是一个可能的值（即Nothing <: T），则应使用Union{Nothing,Some{T}}。例如，findfirst和tryparse使用了这种模式。

这篇博文重点关注第一种情况，希望上面关于missing行为的描述能够清楚地说明为什么将其与nothing区分开来是有用的。事实上，虽然missing在传递给标准数学运算符和函数时通常会传播，但nothing没有实现任何特定方法，因此通常会给出MethodError，迫使调用者显式地处理它。但是，下面讨论的关于性能的考虑同样适用于missing和nothing（以及其他在类似情况下的自定义类型）。

有效的表示

Julia 的另一个优势是，不需要使用诸如向量化调用之类的技巧来使代码快速运行。本着这种精神，处理缺失值必须高效，而无需特殊处理。虽然在以前的 Julia 版本中，Union{Missing,T}方法效率非常低，但由于 Julia 0.7 中编译器实现的两个改进，这种情况发生了巨大变化。

第一个改进涉及对小型Union类型的优化。当类型推断检测到一个变量可以容纳多种类型的值，但这些类型数量有限（例如Union{Missing,T}）时，编译器会为每种可能的类型在单独的分支中生成优化代码，并在检查值的实际类型后运行相应的代码^[6]。这产生的代码非常接近通常与哨兵方法一起使用的代码，在哨兵方法中，需要手动检查处理的值是否等于哨兵。当然，只有当类型被推断为小型Union时，此优化才可用：因此，必须使用Array{Union{Missing,T}}而不是Array{Any}对象，以便为编译器提供必要的类型信息。

第二个改进在于为元素类型是位类型Union的Array对象使用紧凑的内存布局，即不包含引用的不可变类型（参见isbits函数）。这包括Missing和基本类型，如Int、Float64、Complex{Float64}和Date。当T是位类型时，Array{Union{Missing,T}}对象在内部表示为两个相同大小的数组：一个Array{T}保存非缺失值，以及未初始化的缺失值内存；以及一个Array{UInt8}存储类型标签，指示每个条目是Missing类型还是T类型。

这种布局消耗的内存略多于哨兵方法，因为类型标签部分为每个条目占用一个字节。但这种开销是合理的：例如，Array{Union{Missing,Float64}}的内存使用量仅比Array{Float64}高 12.5%。与哨兵方法相比，它具有完全通用的优势（如上所述）。实际上，此机制可用于其他情况，例如Union{Nothing,Int}（这是 Julia 0.7 中indexin返回的数组的元素类型）。

包含缺失值的非位类型数组已经并且继续以与没有缺失值的对应数组一样高效的方式表示。实际上，此类数组包含指向位于不同内存区域的实际对象的指针。缺失值可以像非缺失值一样表示为一个特殊指针。这尤其适用于Array{Union{Missing,String}}。

在存在缺失值的情况下，Array的高效内存布局使得无需使用专用的数组类型，如DataArray。事实上，DataArray类型的布局与上面描述的Array{Union{Missing,T}}非常相似。唯一的区别在于它使用BitArray而不是Array{UInt8}来指示值是否缺失，因此每个条目占用 1 位而不是 8 位。即使它消耗更多的内存，Array{UInt8}掩码方法也更快（至少在当前的BitArray状态下），并且它可以推广到超过两种类型的Union。但是，我们知道其他实现，例如PostgreSQL或Apache Arrow使用等效于BitArray的位图。

这种有效的表示允许编译器在处理包含缺失值的数组时生成非常高效的代码。例如，以下函数是计算数组中非缺失值总和的最直接方法之一

function sum_nonmissing(X::AbstractArray)
    s = zero(eltype(X))
    @inbounds @simd for x in X
        if x !== missing
            s += x
        end
    end
    s
end

在 Julia 0.7 中，当输入为Array{Int}对象时，此相对简单的实现会生成非常高效的本地代码。事实上，由于掩码SIMD指令，缺失值的存在对性能没有任何显著影响。在以下基准测试中，X1是随机Int32条目的向量，X2保存相同的数据，但也可以容纳missing值，X3实际上在随机位置包含 10% 的missing值。sum_nonmissing对所有三个数组花费的时间大致相同（在启用 AVX2 指令的 Intel Skylake CPU 上）。

julia> X1 = rand(Int32, 10_000_000);

julia> X2 = Vector{Union{Missing, Int32}}(X1);

julia> X3 = ifelse.(rand(length(X2)) .< 0.9, X2, missing);

julia> using BenchmarkTools

julia> @btime sum_nonmissing(X1);
  2.738 ms (1 allocation: 16 bytes)

julia> @btime sum_nonmissing(X2);
  3.216 ms (1 allocation: 16 bytes)

julia> @btime sum_nonmissing(X3);
  3.214 ms (1 allocation: 16 bytes)

作为参考点，R 的sum(x, na.rm=TRUE)函数（在 C 中实现）在没有缺失值的情况下大约需要 7 毫秒，在有缺失值的情况下大约需要 19 毫秒（在 R 中，integer数组始终允许缺失值）。

当前限制和未来发展

尽管如此，一切并非完美，仍然需要一些改进。好消息是，最困难的部分已经在 Julia 0.7 中实现了。

第一系列限制涉及性能。即使上面提到的求和示例可以说是令人印象深刻，但在撰写本文时，编译器尚无法在所有情况下生成如此快速的代码。例如，对于Float64元素的数组，缺失值仍然会对性能产生重大影响，而这些数组对于数值计算至关重要

julia> Y1 = rand(10_000_000);

julia> Y2 = Vector{Union{Missing, Float64}}(Y1);

julia> Y3 = ifelse.(rand(length(Y2)) .< 0.9, Y2, missing);

julia> @btime sum_nonmissing(Y1);
  5.733 ms (1 allocation: 16 bytes)

julia> @btime sum_nonmissing(Y2);
  13.854 ms (1 allocation: 16 bytes)

julia> @btime sum_nonmissing(Y3);
  17.780 ms (1 allocation: 16 bytes)

但是，这没有什么可耻的，因为 Julia 仍然比 R 快：对于没有缺失值的numeric数组，sum(x, na.rm=TRUE)大约需要 11 毫秒，对于包含 10% 缺失值的数组大约需要 21 毫秒。这证明了Union{Missing,T}方法的有效性，尽管仍有改进的空间。

编译器改进对于确保即使对于更复杂的模式也能像上面示例中那样生成快速的本地代码也是必要的。特别是，在sum_nonmissing中用ismissing(x)替换x !== missing目前会导致性能大幅下降。还可以注意到，在Array{T}和Array{Union{Missing,T}}之间转换目前涉及复制，理论上可以避免位类型。最后，sum(skipmissing(x))目前由于通用mapreduce实现的工作方式，比自定义sum_nonmissing函数稍微慢一些。希望这些问题在不久的将来能够得到解决，因为它们不需要任何根本性的更改。

另一个可以改进的领域涉及处理缺失值的便捷语法和函数。我们完全意识到，对于那些日常使用缺失值的用户来说，Union{Missing,T}过于冗长。T?语法已被讨论为一种紧凑的替代方案，其灵感来自具有Nullable类型的语言。但是，目前尚不清楚将此语法赋予Union{Missing,T}还是Union{Nothing,T}更合适。因此，它目前被保留，等待做出决定。一个可能的解决方案是为每种类型引入一个专门的语法。

便捷函数对于使用尚未编写为自动传播缺失值的函数传播缺失值也很有用。构造函数如lift(f, x)、lift(f)(x)和f?(x)已被讨论以提供ismissing(x) ? missing : f(x)的更短的等价物。

一个更基本的限制是Union{Missing,T}表示的选择固有的。在此表示中，非缺失值不携带任何关于它是否可能缺失的信息，即关于它是否从允许缺失值的数组或数据集的列中提取的信息。在实践中，这意味着，如果x是Array{Union{Missing,T}}，但实际上不包含缺失值，则map(identity, x)将返回Array{T}。这是因为map仅根据输出的实际内容选择其返回类型，以避免依赖于类型推断（这可能因编译器改进而异）。这也意味着，当将函数应用于Array{Union{T,Missing}}中的每个元素时，无法根据第一个元素的类型选择结果类型，这在例如决定表列是否应允许在 SQL 数据库中使用NULL条目时可能存在问题。这个问题已被多次讨论，但目前尚不清楚哪种缓解方法是最佳的。

尽管存在这些限制，但我们相信 Julia 0.7 中的缺失值支持将是功能和性能方面最完善的，即使是在专门的统计语言中。

作者：Milan Bouchet-Valat，社会学家，法国人口研究所（Ined）研究员，巴黎。

致谢：此框架是多年集体努力的结果。John Myles White 在 2016 年之前领导了围绕 Julia 中缺失值支持的思考。Jameson Nash 和 Keno Fischer 实现了编译器优化，Jacob Quinn 实现了数组的高效内存布局。Alex Arslan、Jeff Bezanson、Stefan Karpinski、Jameson Nash 和 Jacob Quinn 是这项漫长而复杂的设计工作中最主要的参与者。讨论涉及许多其他开发者，有时意见不一，其中 David Anthoff 值得特别提及。

脚注

[1]

C#（使用null）和VB.NET（使用Nothing）是此规则的两个部分例外，因为它们提供了提升的运算符，这些运算符对Nullable参数进行运算并返回Nullable。

[2]

同一软件包中提供的PooledDataArray类型可以用CategoricalArray或PooledArray替换，具体取决于数据是否是真正的分类数据，或者只是包含少量不同的值。

[6]

此优化适用于少量类型的所有Union，无论它们是否为位类型。“少量”由MAX_UNION_SPLITTING常量定义，该常量当前设置为 4。

除了这些promote_rule方法之外，Missing和Nothing类型还实现了内部promote_typejoin函数，该函数确保诸如map和collect之类的函数返回元素类型为Union{Missing,T}或Union{Nothing,T}而不是Any的数组。

[3]

更准确地说，SAS 将缺失值视为小于任何非缺失值，而 Stata 将其视为大于任何非缺失值。这种行为可以通过每种语言的哨兵值的特定选择以及揭示实现细节的不完美抽象来解释。

[4]

例如，参见 Olivier Godechot 的"我们可以使用阿尔萨斯-摩泽尔来估计法国 35 小时工作制法规的就业影响吗？"。

[5]

在2014 年的一篇博文中，John Myles White 提倡使用Nullable，因为它与Union{T,NA}相比具有更高的性能。由于编译器的改进，这种性能差距已不复存在，这些改进彻底改变了缺失值的支持。

Julia 0.7 中一流的统计缺失值支持

2018年6月19日 | Milan Bouchet-Valat

“missing”对象：安全且通用的缺失值

从“Nullable”到“missing”和“nothing”

有效的表示

当前限制和未来发展

脚注