Вызов всех экспертов по локальной регрессии и / или R !
Я столкнулся с ограничением стандартной функции loess в R и надеюсь, что у вас есть совет.Текущая реализация поддерживает только 1-4 предиктора .Позвольте мне изложить наш сценарий применения, чтобы показать, почему это может легко стать проблемой, как только мы захотим использовать глобально подходящие параметрические коваримы .
По существу, у нас есть пространственное искажение s (x, y) с наложением на ряд измерений z :
z_i = s(x_i,y_i) + v_{g_i}
Эти измерения z можно сгруппировать по одному и тому же неискаженному измерениюзначение v для каждой группы g .Членство в группе g_i известно для каждого измерения, но лежащие в основе неискаженные значения измерения v_g для групп неизвестны и должны определяться регрессией ( global , not local).
Нам нужнооценить двумерный пространственный тренд s (x, y) , который мы затем хотим удалить.В нашем приложении, скажем, есть 20 групп по меньшей мере по 35 измерений в наиболее простом сценарии.Измерения размещены случайным образом.Принимая первую группу в качестве эталона, мы получаем 19 неизвестных смещений.
Приведенный ниже код для игрушечных данных (с пространственным трендом в одном измерении x ) работает для двух или трех групп смещений.
К сожалению, вызов loess не выполняется для четырех или более групп смещения с сообщением об ошибке
Error in simpleLoess(y, x, w, span, degree, parametric, drop.square,
normalize, :
only 1-4 predictors are allowed"
Я попытался переопределить ограничение и получил
k>d2MAX in ehg136. Need to recompile with increased dimensions.
Насколько легко это было бы сделать?Я нигде не могу найти определение d2MAX , и кажется, что это может быть жестко закодировано - ошибка, по-видимому, вызывается строкой # 1359 в loessf.f
if(k .gt. 15) call ehg182(105)
В качестве альтернативы кто-нибудь знает о реализации локальной регрессии с глобальными (параметрическими) группами смещений, которая может быть применена здесь?
Или есть лучший способ справиться с этим?Я попробовал lme со структурами корреляции, но это, кажется, намного, намного медленнее.
Любые комментарии будут очень благодарны!
Большое спасибо,
David
###
#
# loess with parametric offsets - toy data demo
#
x<-seq(0,9,.1);
x.N<-length(x);
o<-c(0.4,-0.8,1.2#,-0.2 # works for three but not four
); # these are the (unknown) offsets
o.N<-length(o);
f<-sapply(seq(o.N),
function(n){
ifelse((seq(x.N)<= n *x.N/(o.N+1) &
seq(x.N)> (n-1)*x.N/(o.N+1)),
1,0);
});
f<-f[sample(NROW(f)),];
y<-sin(x)+rnorm(length(x),0,.1)+f%*%o;
s.fs<-sapply(seq(NCOL(f)),function(i){paste('f',i,sep='')});
s<-paste(c('y~x',s.fs),collapse='+');
d<-data.frame(x,y,f)
names(d)<-c('x','y',s.fs);
l<-loess(formula(s),parametric=s.fs,drop.square=s.fs,normalize=F,data=d,
span=0.4);
yp<-predict(l,newdata=d);
plot(x,y,pch='+',ylim=c(-3,3),col='red'); # input data
points(x,yp,pch='o',col='blue'); # fit of that
d0<-d; d0$f1<-d0$f2<-d0$f3<-0;
yp0<-predict(l,newdata=d0);
points(x,y-f%*%o); # spatial distortion
lines(x,yp0,pch='+'); # estimate of that
op<-sapply(seq(NCOL(f)),function(i){(yp-yp0)[!!f[,i]][1]});
cat("Demo offsets:",o,"\n");
cat("Estimated offsets:",format(op,digits=1),"\n");