Быстрый взвешенный расчет матрицы рассеяния

В этом вопросе шесть месяцев назад , jez был достаточно хорош, чтобы помочь мне придумать быстрое приближение для внешнего произведения разностей строк, то есть:

K = np.zeros((len(X), len(X))) for i, Xi in enumerate(X): for j, Xj in enumerate(X): dij = Xi - Xj K += np.outer(dij, dij) 

Это помогло найти расчет матрицы рассеяния для формы дискриминантного анализа Фишера. Но теперь я пытаюсь сделать локальный анализ дискриминантов Fisher, где каждый внешний продукт взвешен матрицей A, которая имеет информацию о местонахождении пары, поэтому новая строка:

K += A[i][j] * np.outer(dij, dij)

К сожалению, быстрый способ рассчитать невзвешенную матрицу рассеяния, представленную в предыдущем ответе, не работает для этого, и, насколько я могу судить, сделать быстрые изменения непросто.

Линейная алгебра, безусловно, не мой сильный костюм, я не умею придумывать подобные вещи. Что такое быстрый способ вычисления взвешенной суммы парных внешних разностных продуктов?

    One Solution collect form web for “Быстрый взвешенный расчет матрицы рассеяния”

    Вот способ векторизации указанного вами расчета. Если вы так много делаете, то, возможно, стоит изучить, как использовать «numpy.tensordot». Он умножает все элементы в соответствии со стандартным трансляцией numpy, а затем суммирует по парам осей, заданных с помощью kwrd, «axes».

    Вот код:

     # Imports import numpy as np from numpy.random import random # Original calculation for testing purposes def ftrue(A, X): "" K = np.zeros((len(X), len(X))) KA_true = np.zeros((len(X), len(X))) for i, Xi in enumerate(X): for j, Xj in enumerate(X): dij = Xi - Xj K += np.outer(dij, dij) KA_true += A[i, j] * np.outer(dij, dij) return ftrue # Better: No Python loops. But, makes a large temporary array. def fbetter(A, X): "" c = X[:, None, :] - X[None, :, :] b = A[:, :, None] * c # ! BAD ! temporary array size N**3 KA_better = np.tensordot(b, c, axes = [(0,1),(0,1)]) return KA_better # Best way: No Python for loops. No large temporary arrays def fbest(A, X): "" KA_best = np.tensordot(A.sum(1)[:,None] * X, X, axes=[(0,), (0,)]) KA_best += np.tensordot(A.sum(0)[:,None] * X, X, axes=[(0,), (0,)]) KA_best -= np.tensordot(np.dot(A, X), X, axes=[(0,), (0,)]) KA_best -= np.tensordot(X, np.dot(A, X), axes=[(0,), (0,)]) return KA_best # Test script if __name__ == "__main__": # Parameters for the computation N = 250 X = random((N, N)) A = random((N, N)) # Print the error KA_better = fbetter(A, X) KA_best = fbest(A, X) # Test against true if array size isn't too big if N<100: KA_true = ftrue(A, X) err = abs(KA_better - KA_true).mean() msg = "Mean absolute difference (better): {}." print(msg.format(err)) # Test best against better err = abs(KA_best - KA_better).mean() msg = "Mean absolute difference (best): {}." print(msg.format(err)) 

    Моя первая попытка (fbetter) сделала большой временный массив размера NxNxN. Вторая попытка (fbest) никогда не делает ничего большего, чем NxN. Это работает довольно хорошо до N ~ 1000.

    Тест времени

    Кроме того, код работает быстрее, когда выходной массив меньше. введите описание изображения здесь

    У меня установлен MKL, поэтому призывы к tensordot действительно быстрые и работают параллельно.

    Спасибо за вопрос. Это было отличное упражнение и напомнило мне, насколько важно избегать создания больших временных массивов.

     
    Interesting Posts for Van-Lav

    Python – Как открыть файл и указать смещение в байтах?

    Параллельность с модулем подпроцесса. Как я могу это сделать?

    Запуск команды Судо с парамико

    Ошибка при запуске носетистов

    Повторная выборка с «how = count» вызывает проблемы

    Использование gevent monkey patching с резьбой делает работу нити серийно

    Буферные объявления Cython для участников объекта

    Что такое кросс-платформенный способ воспроизведения звукового файла в python?

    Ошибка Python / Flask: «ImportError: невозможно импортировать имя _compare_digest»

    Преобразование строки в значение ASCII python

    Развертывание matplotlib на heroku не удалось. Как это сделать правильно?

    Как обновить метаданные существующего объекта в AWS S3 с помощью python boto3?

    Создание автоматических тестов для интерактивной оболочки на основе модуля cmd на Python

    Создание и поддержка нескольких сеансов ssh в Python

    двоичная логистическая регрессия xgboost

    Python - лучший язык программирования в мире.