Как бы вы сделали эту матричную операцию в MATLAB?

Question

Учитывая две случайные величины X и Y, где X = (x1, .., xn) и Y = (y1, ..., yn) в nx2 матрице A, поэтому A=[X Y], мне нужно выполните следующую операцию:

медиана ((x-медиана (x)) (y-медиана (y)))

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

Пока Я сделал следующий код:

 for i=1:n
      for j=1:n
         a1=median(A(:,i));
         a2=median(A(:,j));
         SMM(i,j)=median(((A(:,i)-a1(ones(t,1),:)).*(A(:,j)-a2(ones(t,1),:))));
      end
 end 

Однако теоретически я должен получить полуопределенную (положительную или отрицательную) симметричную c матрицу, однако в этом коде дело обстоит иначе.

Я ошибаюсь в формулировке кода?

Answer 1

Различные точки:

• Для каждого из ваших столбцов A ( x , y ), медиана (a1, a2) не меняется. Вы должны вычислить их вне циклов.

• Петли go более n, а не t, которые являются переменными и индексами выходной матрицы.

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

A = A - median(A,1);  % be explicit about which dimension to take the median over!

Далее, мы бы l oop над t x t выводит элементы ковариационной матрицы и вычисляет каждый из элементов:

t = size(A,2);
SMM = zeros(t,t);  % always preallocate output arrays before a loop
for j=1:t
  for i=1:t
    SMM(i,j) = median(A(:,i).*A(:,j));
  end
end

l oop может быть векторизовано, но это приводит к большой промежуточной матрице, которая замедляется код также. Так что это может не стоить усилий, чтобы векторизовать. Попробуйте, только если этот код слишком медленный!

Также должна быть возможность запустить внутреннюю l oop из i=j:t, чтобы пропустить вычисление избыточной половины матрицы симметрии c вместо копирования по ранее вычисленным значениям.