Чрезвычайно большое средневзвешенное

Question

Я использую 64-битный Matlab с 32 г оперативной памяти (как вы знаете).

У меня есть файл (вектор) из 1,3 миллиона чисел (целых чисел). Я хочу сделать еще один вектор такой же длины, где каждая точка представляет собой средневзвешенное значение всего первого вектора, взвешенное на обратном расстоянии от этой позиции (на самом деле это позиция ^ -0,1, а не ^ -1, но для примера) , Я не могу использовать функцию фильтра 'matlab', потому что она может только усреднять вещи до текущей точки, верно? Чтобы объяснить более четко, вот пример из 3 элементов

data = [ 2 6 9 ]
weights = [ 1 1/2 1/3; 1/2 1 1/2; 1/3 1/2 1 ]
results=data*weights= [ 8 11.5 12.666 ]
i.e.
8 = 2*1 + 6*1/2 + 9*1/3
11.5 = 2*1/2 + 6*1 + 9*1/2
12.666 = 2*1/3 + 6*1/2 + 9*1

Таким образом, каждая точка в новом векторе является средневзвешенным значением всего первого вектора, взвешиваясь на 1 / (расстояние от этой позиции + 1).

Я мог бы просто переделать вектор весов для каждой точки, а затем вычислить вектор результатов элемент за элементом, но для этого требуется 1,3 миллиона итераций цикла for, каждая из которых содержит 1,3 миллиона умножений. Я бы предпочел использовать прямое матричное умножение, умножив 1x1.3mil на 1.3milx1.3mil, что теоретически работает, но я не могу загрузить такую ​​большую матрицу.

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

Мне не нужно делать это в matlab, поэтому любые советы, которые люди дают по поводу использования таких больших чисел и получения средних значений, были бы признательны. Поскольку я использую вес ^ -0,1, а не ^ -1, он не падает так быстро - миллионная точка по-прежнему весит 0,25 по сравнению с исходным весом точек 1, поэтому я не могу просто сократить его и он становится большим.

Надеюсь, это было достаточно ясно?

Вот код ответа ниже (чтобы его можно было отформатировать?):

data = load('/Users/mmanary/Documents/test/insertion.txt');
data=data.';
total=length(data);
x=1:total;
datapad=[zeros(1,total) data];
weights = ([(total+1):-1:2 1:total]).^(-.4);
weights = weights/sum(weights);
Fdata = fft(datapad);
Fweights = fft(weights);
Fresults = Fdata .* Fweights;
results = ifft(Fresults);
results = results(1:total);
plot(x,results)

Answer 1

Единственный разумный способ сделать это - использовать FFT сверточный , так как он лежит в основе функции filter и аналогичных.Это очень легко сделать вручную:

% Simulate some data
n = 10^6;
x = randi(10,1,n);
xpad = [zeros(1,n) x];

% Setup smoothing kernel
k = 1 ./ [(n+1):-1:2 1:n];

% FFT convolution
Fx = fft(xpad);
Fk = fft(k);

Fxk = Fx .* Fk;

xk = ifft(Fxk);
xk = xk(1:n);

Занимает менее чем за полсекунды для n=10^6!

Answer 2

Это, вероятно, не лучший способ сделать это, но с большим объемом памяти вы можете определенно распараллелить процесс.

Вы можете построить разреженные матрицы, состоящие из элементов вашей исходной матрицы, которые имеют значение i^(-1) (где i = 1 .. 1.3 million), умножить их на исходный вектор и суммировать все результаты вместе.

Так что для вашего примера продукт будет по существу:

a = rand(3,1);
b1 = [1 0 0;
      0 1 0;
      0 0 1];
b2 = [0 1 0;
      1 0 1;
      0 1 0] / 2;
b3 = [0 0 1;
      0 0 0;
      1 0 0] / 3;

c = sparse(b1) * a + sparse(b2) * a + sparse(b3) * a;

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

Посмотрите на цикл parfor в MATLAB: http://www.mathworks.com/help/toolbox/distcomp/parfor.html

Answer 3

Я не могу использовать функцию фильтра 'matlab', потому что она может только усреднять вещи до текущей точки, верно?

Это не правильно.Вы всегда можете добавить образцы (т.е. добавить или удалить нули) из ваших данных или из отфильтрованных данных.Поскольку фильтрация с помощью filter (кстати, вы также можете использовать conv) - это линейное действие, оно не изменит результат (это как добавление и удаление нулей, которые ничего не делают, а затем фильтрация. Тогда линейность позволяет вампоменять местами порядок добавления образцов -> фильтр -> удалить образец).

В любом случае, в вашем примере вы можете взять ядро ​​усреднения следующим образом:

weights = 1 ./ [3 2 1 2 3]; % this kernel introduces a delay of 2 samples

, а затем просто:

result =  filter(w,1,[data, zeros(1,3)]); % or conv (data, w)
% removing the delay introduced by the kernel
result = result (3:end-1);

Answer 4

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

Операции ^ -0.1 для создания весов будут намного дольше, чем операции + и * для вычислениясредневзвешенное значение, но вы повторно используете весовые коэффициенты для всех миллионов средневзвешенных операций.Алгоритм становится:

• Создать вектор весов со всеми весами, необходимыми для любого вычисления: weights = (-n:n).^-0.1

• Для каждого элемента в векторе:

  • Индексируйте соответствующую часть вектора weights, чтобы считать текущий элемент «центром».

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

Основной цикл выполняет n ^ 2 сложения иsubractions.При n , равном 1,3 млн, это 3,4 трлн операций.Одно ядро ​​современного процессора 3GHz может делать 6 миллиардов сложений / умножений в секунду, так что получается около 10 минут.Добавьте время для индексации вектора weights и накладных расходов, и я все еще рассчитываю, что вы можете прийти через полчаса.

Answer 5

Вы рассматривали только 2 варианта: Умножение 1,3М * 1,3М матрицы на вектор один раз или умножение 2 1,3М векторов в 1,3М раз.

Но вы можете разделить свою матрицу весов на столько подматриц, сколько пожелаете, и умножить матрицу n * 1,3M на вектор 1,3M / n раз.

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

с вашим объемом памяти вы должны начать с n = 5000.

вы также можете сделать это быстрее, используя parfor (с n, деленным на количество процессоров).