Регрессия в R с использованием векторизации и матриц - PullRequest
0 голосов
/ 23 января 2012

У меня есть векторизация Q в R с использованием матриц. У меня есть 2 Cols, которые нужно регрессировать против каждого, используя определенные индексы. Данные

matrix_senttoR = [ ...
                  0.11 0.95
                  0.23 0.34
                  0.67 0.54
                  0.65 0.95
                  0.12 0.54
                  0.45 0.43 ] ;
indices_forR = [ ...
            1
            1
            1
            2
            2
            2 ] ;

Col1 в матрице - это данные, скажем, MSFT и GOOG (по 3 строки в каждой), а Col2 - доход от эталонного StkIndex на соответствующие даты. Данные в матричном формате, так как они отправлены из Matlab.

Я сейчас использую

slope <- by(    data.frame(matrix_senttoR),   indices_forR,   FUN=function(x)  
                         {zyp.sen (X1~X2,data=x) $coeff[2] }      ) 
betasFac <- sapply(slope , function(x) x+0)

Я использую data.frame выше, так как я не мог использовать cbind (). Если я использую cbind (), то Matlab выдает ошибку, поскольку не понимает этот формат данных. Я запускаю эти команды изнутри Matlab (http://www.mathworks.com/matlabcentral/fileexchange/5051). Вы можете заменить zyp (zyp.sen) на lm.

BY здесь медленно (может быть из-за фреймов данных?). Есть ли лучший способ сделать это? Требуется 14 сек + для 150 тыс. Строк данных. Могу ли я вместо этого использовать матричную векторизацию в R? Спасибо.

Ответы [ 2 ]

1 голос
/ 23 января 2012

Это можно легко переместить в комментарий, но:

Несколько вещей, на которые стоит обратить внимание, я стараюсь избегать функции by(), так как ее возвращаемое значение является фанки-объектом. Вместо этого, почему бы не добавить свой вектор indices_forR в data.frame?

df <- data.frame(matrix_senttoR) 
df$indices_forR <- indices_forR

пакет plyr выполняет работу отсюда:

ddply(df,.(indices_forR),function(x) zyp.sen(X1~X2,data=x)$coeff[2])

Вы можете легко многопоточность эту операцию, используя doMC или doSnow и аргумент .parallel=TRUE для ddply.

если цель - скорость, я бы также изучил пакет data.table (который упаковывает data.frame и работает намного быстрее). Кроме того, я предполагаю, что медленная часть - это вызов zyp.sen(), а не by(). Выполнение на нескольких ядрах ускорит это.

> dput(df)
structure(list(X1 = c(0.11, 0.23, 0.67, 0.65, 0.12, 0.45), X2 = c(0.95, 
0.34, 0.54, 0.95, 0.54, 0.43), indices_forR = c(1, 1, 1, 2, 2, 
2)), .Names = c("X1", "X2", "indices_forR"), row.names = c(NA, 
-6L), class = "data.frame")

> ddply(df,.(indices),function(x) lm(X1~X2,data=x)$coeff[2])
  indices         X2
1       1 -0.3702172
2       2  0.6324900
0 голосов
/ 24 января 2012

Я все еще думаю, что вы слишком усложняете вещи, переходя от MATLAB к R и обратно. А передача 150 тыс. Строк данных должна значительно замедлить процесс.

zyp.sen на самом деле довольно тривиально портировать на MATLAB. Вот, пожалуйста:

function [intercept, slope, intercepts, slopes, rank, residuals] = ZypSen(x, y)
% Computes a Thiel-Sen estimate of slope for a vector of data.

n = length(x);

slopes = arrayfun(@(i) ZypSlopediff(i, x, y, n), 1:(n - 1), ...
    'UniformOutput', false);
slopes =  [slopes{:}];
sni = isfinite(slopes);
slope = median(slopes(sni));

intercepts = y - slope * x;
intercept = median(intercepts);

rank = 2;
residuals = x - slope * y + intercept;

end


function z = ZypSlopediff(i, x, y, n)

z = (y(1:(n - i)) - y((i + 1):n)) ./ ...
    (x(1:(n - i)) - x((i + 1):n));

end

Я проверил это, используя R example(zyp.sen), и он дает тот же ответ.

x = [0 1 2 4 5]
y = [6 4 1 8 7]
[int, sl, ints, sls, ra, res] = ZypSen(x, y)

Вы действительно должны провести дополнительную проверку, просто чтобы быть уверенным.

...