ROCKET-MiniROCKET-MultiROCKET代码解析

源码

ROCKET: exceptionally fast and accurate time series classification using random convolutional kernels

MiniRocket: A Very Fast (Almost) Deterministic Transform for Time Series Classification

MultiRocket: multiple pooling operators and
transformations for fast and effective time series
classification

ROCKET

ROCKET在sktime包种可直接调用
ROCKET使用教程
内容总结：
使用随机卷积核作为时间序列特征提取方法，使用随机卷积核长度、dilation、padding，随机初始化，不需要学习权重;
stride=1；
不使用任何非线性（激活函数），使用ppv作为全局最大池化方法；
偏置作为ppv的阈值；
对于ROCKET，只需修改卷积核个数，一般来说，K越大精确度越高但是时间也越长；
每个卷积核对时间序列进行应用得到特征图，ROCKET 从每个特征图中计算两个值（最大值、特征图中正值的比例）作为输出的特征值。

generate_kernels

#随机卷积核参数
def generate_kernels(input_length, num_kernels):

    candidate_lengths = np.array((7, 9, 11), dtype = np.int32)#卷积核长度范围7，9，11
    lengths = np.random.choice(candidate_lengths, num_kernels)#对每个卷积核随机选择卷积核长度
    #创建weight，biase，dilation,padding数组
    weights = np.zeros(lengths.sum(), dtype = np.float64)
    biases = np.zeros(num_kernels, dtype = np.float64)
    dilations = np.zeros(num_kernels, dtype = np.int32)
    paddings = np.zeros(num_kernels, dtype = np.int32)

    a1 = 0

    for i in range(num_kernels):

        _length = lengths[i]

        _weights = np.random.normal(0, 1, _length)#权重初始化参数服从正态分布

        b1 = a1 + _length
        weights[a1:b1] = _weights - _weights.mean()

        biases[i] = np.random.uniform(-1, 1)#偏执初始化参数服从均匀分布

        dilation = 2 ** np.random.uniform(0, np.log2((input_length - 1) / (_length - 1)))#dilation初始化服从2的指数
        dilation = np.int32(dilation)
        dilations[i] = dilation

        padding = ((_length - 1) * dilation) // 2 if np.random.randint(2) == 1 else 0#根据随机数决定是否填充
        paddings[i] = padding

        a1 = b1

    return weights, lengths, biases, dilations, paddings

apply_kernel

计算一个卷积核输出的特征值

def apply_kernel(X, weights, length, bias, dilation, padding):

    input_length = len(X)

    output_length = (input_length + (2 * padding)) - ((length - 1) * dilation)#输出的特征图长度

    _ppv = 0
    _max = np.NINF

    end = (input_length + padding) - ((length - 1) * dilation)

    for i in range(-padding, end):

        _sum = bias

        index = i

        for j in range(length):

            if index > -1 and index < input_length:

                _sum = _sum + weights[j] * X[index]#计算当前的特征值

            index = index + dilation

        if _sum > _max:
            _max = _sum#计算最大值

        if _sum > 0:
            _ppv += 1#计算正值个数

    return _ppv / output_length, _max

apply_kernels

计算所有卷积特征值

def apply_kernels(X, kernels):

    weights, lengths, biases, dilations, paddings = kernels

    num_examples, _ = X.shape
    num_kernels = len(lengths)

    _X = np.zeros((num_examples, num_kernels * 2), dtype = np.float64) # 创建输出特征数组，每个卷积核有两个特征值

    for i in prange(num_examples):

        a1 = 0 # for weights
        a2 = 0 # for features

        for j in range(num_kernels):

            b1 = a1 + lengths[j]
            b2 = a2 + 2

            _X[i, a2:b2] = \
            apply_kernel(X[i], weights[a1:b1], lengths[j], biases[j], dilations[j], paddings[j])#记录特征值

            a1 = b1#更新索引
            a2 = b2

    return _X

MiniROCKET

MiniROCKET在sktime包种可直接调用
MiniROCKET使用示例
内容总结：
Mini ROCKET提高了计算速度；
将卷积核长度固定为9;
将卷积核权重限制在两个值（-1，2），偏置值是从卷积输出中提取的，对于一个随机选择的训练样本，给定一个卷积核，计算其输出 ，然后取[0.25, 0.5, 0.75]分位数作为偏置的可选值;
dilation的值确定范围为[2⁰, 2^max]，max=log2(l(input)-1)/(l(kernel length)-1);
Padding: 一半使用填充，一半不使用填充，使用zero padding
特征：仅使用PPV，不考虑最大值
特征值个数=卷积核个数(84)×dilations数×bias数
之所以是84个卷积核是因为卷积核个数为9，其中6个权重为-1，3个权重为2，它的所有排列组合是84个组合

_fit_dilations

计算一个卷积核需要的dilation数以及一个dilation生成的特征值个数

def _fit_dilations(n_timepoints, num_features, max_dilations_per_kernel):
    num_kernels = 84#卷积核子集个数

    num_features_per_kernel = num_features // num_kernels#每个卷积核需要提供的特征个数
    true_max_dilations_per_kernel = min(
        num_features_per_kernel, max_dilations_per_kernel
    )#每个卷积核的dalition个数
    multiplier = num_features_per_kernel / true_max_dilations_per_kernel#每个dalition给每个卷积核提供的特征数

    max_exponent = np.log2((n_timepoints - 1) / (9 - 1))#最大膨胀指数
    dilations, num_features_per_dilation = np.unique(
        np.logspace(0, max_exponent, true_max_dilations_per_kernel, base=2).astype(
            np.int32
        ),
        return_counts=True,
    )#随机生成dalition值
    num_features_per_dilation = (num_features_per_dilation * multiplier).astype(
        np.int32
    )  #每个dalition给每个卷积核提供的特征数
    remainder = num_features_per_kernel - np.sum(num_features_per_dilation)
    i = 0
    while remainder > 0:
        num_features_per_dilation[i] += 1
        remainder -= 1
        i = (i + 1) % len(num_features_per_dilation)

    return dilations, num_features_per_dilation

_fit_biases_multi

生成bias值

def _fit_biases_multi(
    X,
    num_channels_per_combination,
    channel_indices,
    dilations,
    num_features_per_dilation,
    quantiles,
    seed,
):
    if seed is not None:
        np.random.seed(seed)

    n_instances, n_columns, n_timepoints = X.shape

    # equivalent to:
    from itertools import combinations
    indices = np.array([_ for _ in combinations(np.arange(9), 3)]) #卷积核β权重的位置，每个卷积核有6个-1，3个2

    num_kernels = len(indices)
    num_dilations = len(dilations)

    num_features = num_kernels * np.sum(num_features_per_dilation)

    biases = np.zeros(num_features, dtype=np.float32)

    feature_index_start = 0

    combination_index = 0
    num_channels_start = 0

    for dilation_index in range(num_dilations):
        dilation = dilations[dilation_index]
        padding = ((9 - 1) * dilation) // 2

        num_features_this_dilation = num_features_per_dilation[dilation_index]#这个膨胀因子需要的特征值

        for kernel_index in range(num_kernels):
            feature_index_end = feature_index_start + num_features_this_dilation#特征值对应的索引

            num_channels_this_combination = num_channels_per_combination[
                combination_index
            ]#这个卷积核的有效通道数

            num_channels_end = num_channels_start + num_channels_this_combination

            channels_this_combination = channel_indices[
                num_channels_start:num_channels_end
            ]#对应的有效通道

            _X = X[np.random.randint(n_instances)][channels_this_combination]#随机选择一组样本保留有效通道数据

            A = -_X #用于卷积核权重计算 # A = alpha * X = -X
            G = _X + _X + _X  # G = gamma * X = 3X

            C_alpha = np.zeros(
                (num_channels_this_combination, n_timepoints), dtype=np.float32
            )
            C_alpha[:] = A

            C_gamma = np.zeros(
                (9, num_channels_this_combination, n_timepoints), dtype=np.float32
            )
            C_gamma[9 // 2] = G#向下取整除

            start = dilation
            end = n_timepoints - padding

            for gamma_index in range(9 // 2):
                C_alpha[:, -end:] = C_alpha[:, -end:] + A[:, :end]
                C_gamma[gamma_index, :, -end:] = G[:, :end]

                end += dilation

            for gamma_index in range(9 // 2 + 1, 9):
                C_alpha[:, :-start] = C_alpha[:, :-start] + A[:, start:]
                C_gamma[gamma_index, :, :-start] = G[:, start:]

                start += dilation

            index_0, index_1, index_2 = indices[kernel_index]

            C = C_alpha + C_gamma[index_0] + C_gamma[index_1] + C_gamma[index_2]#把卷积乘转换成了加法
            C = np.sum(C, axis=0)

            biases[feature_index_start:feature_index_end] = np.quantile(
                C, quantiles[feature_index_start:feature_index_end]
            )#计算bias

            feature_index_start = feature_index_end

            combination_index += 1
            num_channels_start = num_channels_end

    return biases

_transform_multi

求特征值，特征值个数为84的整数倍，只能接近用户设置的特征值个数

def _transform_multi(X, parameters):
    n_instances, n_columns, n_timepoints = X.shape

    (
        num_channels_per_combination,
        channel_indices,
        dilations,
        num_features_per_dilation,
        biases,
    ) = parameters

    # equivalent to:

    indices = np.array([_ for _ in combinations(np.arange(9), 3)])


    num_kernels = len(indices)
    num_dilations = len(dilations)

    num_features = num_kernels * np.sum(num_features_per_dilation)

    features = np.zeros((n_instances, num_features), dtype=np.float32)
    #计算特征值
    for example_index in prange(n_instances):#变换样本
        _X = X[example_index]

        A = -_X  # A = alpha * X = -X
        G = _X + _X + _X  # G = gamma * X = 3X

        feature_index_start = 0

        combination_index = 0
        num_channels_start = 0

        for dilation_index in range(num_dilations):#变换膨胀因子
            _padding0 = dilation_index % 2

            dilation = dilations[dilation_index]
            padding = ((9 - 1) * dilation) // 2

            num_features_this_dilation = num_features_per_dilation[dilation_index]

            C_alpha = np.zeros((n_columns, n_timepoints), dtype=np.float32)
            C_alpha[:] = A

            C_gamma = np.zeros((9, n_columns, n_timepoints), dtype=np.float32)
            C_gamma[9 // 2] = G

            start = dilation
            end = n_timepoints - padding
            #根据膨胀值来构建矩阵
            for gamma_index in range(9 // 2):
                C_alpha[:, -end:] = C_alpha[:, -end:] + A[:, :end]
                C_gamma[gamma_index, :, -end:] = G[:, :end]

                end += dilation

            for gamma_index in range(9 // 2 + 1, 9):
                C_alpha[:, :-start] = C_alpha[:, :-start] + A[:, start:]
                C_gamma[gamma_index, :, :-start] = G[:, start:]

                start += dilation

            for kernel_index in range(num_kernels):#变换卷积核权重
                feature_index_end = feature_index_start + num_features_this_dilation

                num_channels_this_combination = num_channels_per_combination[
                    combination_index
                ]

                num_channels_end = num_channels_start + num_channels_this_combination

                channels_this_combination = channel_indices[
                    num_channels_start:num_channels_end
                ]#有效通道

                _padding1 = (_padding0 + kernel_index) % 2

                index_0, index_1, index_2 = indices[kernel_index]

                C = (
                    C_alpha[channels_this_combination]
                    + C_gamma[index_0][channels_this_combination]
                    + C_gamma[index_1][channels_this_combination]
                    + C_gamma[index_2][channels_this_combination]
                )
                C = np.sum(C, axis=0)#卷积计算
                #根据不同的bias计算特征值，因为bias可以看成是一个阈值，阈值不同，ppv不同
                if _padding1 == 0:
                    for feature_count in range(num_features_this_dilation):
                        features[
                            example_index, feature_index_start + feature_count
                        ] = _PPV(C, biases[feature_index_start + feature_count]).mean()
                else:
                    for feature_count in range(num_features_this_dilation):
                        features[
                            example_index, feature_index_start + feature_count
                        ] = _PPV(
                            C[padding:-padding],
                            biases[feature_index_start + feature_count],
                        ).mean()

                feature_index_start = feature_index_end

                combination_index += 1
                num_channels_start = num_channels_end

    return features

MultiRocket

主要内容：
使用与 MiniRocket 相同的内核集；
MultiRocket 增加了原始数据的一阶差分，由于长度不同，一阶差分与原始数据拥有不同的偏置和dilation；
增加了 3 个额外的池化算子，正值均值 （MPV）、正值位置均值 （MIPV） 和最长正值延伸 （LSPV）；

transform

基本和MiniRocket差不多，就是再多加了正值均值 （MPV）、正值位置均值 （MIPV） 和最长正值延伸 （LSPV）这几个特征值的计算以及差分数据的特征值计算

                    if _padding1 == 0:
                        for feature_count in range(num_features_this_dilation):
                            feature_index = feature_index_start + feature_count
                            _bias = biases[feature_index]

                            ppv = 0#正值个数
                            last_val = 0#上一次负值位置
                            max_stretch = 0.0#最长连续正值长度
                            mean_index = 0#正值位置累加
                            mean = 0#正值大小累加

                            for j in range(C.shape[0]):
                                if C[j] > _bias:
                                    ppv += 1
                                    mean_index += j
                                    mean += C[j] + _bias
                                elif C[j] < _bias:
                                    stretch = j - last_val
                                    if stretch > max_stretch:
                                        max_stretch = stretch
                                    last_val = j
                            stretch = C.shape[0] - 1 - last_val
                            if stretch > max_stretch:
                                max_stretch = stretch

                            end = feature_index
                            features[example_index, end] = ppv / C.shape[0]
                            end = end + num_features
                            features[example_index, end] = max_stretch
                            end = end + num_features
                            features[example_index, end] = mean / ppv if ppv > 0 else 0
                            end = end + num_features
                            features[example_index, end] = mean_index / ppv if ppv > 0 else -1
                    else:
                        _c = C[padding:-padding]

                        for feature_count in range(num_features_this_dilation):#一个dilation中的特征
                            feature_index = feature_index_start + feature_count
                            _bias = biases[feature_index]

                            ppv = 0#正值个数
                            last_val = 0#上一次负值位置
                            max_stretch = 0.0#最长连续正值长度
                            mean_index = 0#正值位置累加
                            mean = 0#正值大小累加

                            for j in range(_c.shape[0]):
                                if _c[j] > _bias:
                                    ppv += 1
                                    mean_index += j
                                    mean += _c[j] + _bias
                                elif _c[j] < _bias:
                                    stretch = j - last_val
                                    if stretch > max_stretch:
                                        max_stretch = stretch
                                    last_val = j
                            stretch = _c.shape[0] - 1 - last_val
                            if stretch > max_stretch:
                                max_stretch = stretch

                            end = feature_index
                            features[example_index, end] = ppv / _c.shape[0]
                            end = end + num_features
                            features[example_index, end] = max_stretch
                            end = end + num_features
                            features[example_index, end] = mean / ppv if ppv > 0 else 0
                            end = end + num_features
                            features[example_index, end] = mean_index / ppv if ppv > 0 else -1