R Сравнение всех столбцов в матрице от каждого в цикле

Question

У меня есть матрица переменных, и я пытаюсь запустить цикл, который сравнивает разницу между всеми переменными в регрессии, так что создается и заполняется матрица разница. Ниже приведен некоторый код моделирования для решения проблемы. Я хотел бы создать матрицу, которая сравнивает x_1, x_2 и x_3, чтобы создать матрицу 3x3, которая симметрична относительно диагонали, и все они должны быть нулями.

y <- sample(seq(1:4), 100, replace = TRUE) 
x_1 <- sample(seq(1:2), 100, replace = TRUE) 
x_2 <- sample(seq(1:4), 100, replace = TRUE) 
x_3 <- sample(seq(1:4), 100, replace = TRUE) 

frame <- cbind(x_1, x_2, x_3) 
dif <- matrix(NA, ncol = 3, nrow = 3) 

for(i in 1:3){ 
    model_1 <- lm(y ~ frame[,i]) 
    model_2 <- lm(y ~ frame[,i]) 
    dif[i]<- (model_2$coef[2] - model_1$coef[2]) 
} 

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

Answer 1

vcoef <- numeric(3) 
for(i in 1:3) { 
    vcoef[i] <- coef(lm(y~frame[,i]))[2] 
       } 

outer(vcoef, vcoef, "-") 
#---------- 
      [,1]  [,2]  [,3] 
[1,] 0.0000000 -0.15208933 -0.17302592 
[2,] 0.1520893 0.00000000 -0.02093659 
[3,] 0.1730259 0.02093659 0.00000000 

Если вы не хотите избыточную информацию, которую вы могли бы получить все попарные различия с combn:

> combcos <- combn(vcoef,2) 
> combcos[1, ] -combcos[2, ] 
[1] -0.15208933 -0.17302592 -0.02093659 

Answer 2

Попробуйте это:

model <- list() 
for(i in 1:3) { 
    model[[i]] <- lm(y~frame[,i]) 
} 

dif<-sapply(1:3, function(i) { sapply(1:3, function(j) { model[[i]]$coef[2] - model[[j]]$coef[2] }) }) 

матрица будет антисимметричен, т.е. dif[i,j] = -dif[j,i]

Answer 3

Я не могу понять, что вы после этого, но вот трещины с использованием версии outer, которая принимает множество векторов.

library(qdap) 
FUN <- function(x1, x2, y)lm(y ~ x1)$coef[2] - lm(y ~ x2)$coef[2] 
v_outer(list(x_1, x_2, x_3), FUN, y = y) 

##  X1  X2  X3 
## X1 0.000 -0.311 -0.079 
## X2 0.311 0.000 0.232 
## X3 0.079 -0.232 0.000 

Answer 4

Я предпочитаю маршрут eval и parse и как @Tyler мне нравится base:::outer ...

# Make your data into a data.frame 
df <- data.frame(y , x_1 , x_2 , x_3) 

# The variables we want to test 
x <- c("x_1","x_2","x_3") 

# Make the text for each model to parse and evalaute 
mods <- paste0("lm(y ~ " , x , " , data = df)") 

# Evaluate the lm for each variable 
coefs <- unlist(lapply(mods , function(x) eval(parse(text=x))$coef[2])) 
#  x_1   x_2   x_3 
# -0.52140856 0.04662379 0.08694344 

# Combine the results with outer 
outer(coefs , coefs , "-") 
#   x_1   x_2   x_3 
# x_1 0.0000000 -0.56803236 -0.60835201 
# x_2 0.5680324 0.00000000 -0.04031965 
# x_3 0.6083520 0.04031965 0.00000000 

Answer 5

Если я понимаю вы пытаетесь сравнить модели, сравнивая их коэффициенты. Одна из идей - использовать пакет meifly.

Сначала я сгенерировать данные:

set.seed(1) 
frame <- matrix(sample(1:4,3*100,rep=TRUE),ncol=3) 
y <- sample(seq(1:4), 100, replace = TRUE) 

Затем я использую fitbest, который использует leaps пакет очень быстро найти п лучшие модели для заданного числа переменных.

library(meifly) 
library(reshape2) 
library(ggplot2) 
## we look only on models with one variable 
res <- fitbest(y~.,as.data.frame(frame),nvmax=1) 
## get coefficients 
res.coef <- coef(res) 
## remove zero models 
res.coef[res.coef == 0] <- NA 
res.coef <- na.omit(res.coef) 

Теперь для каждой модели имеется сводка коэффициентов. Для каждой переменной мы имеем следующую информацию:

1. Сырое коэффициент
2. т значение и абсолютное значение т-
3. стандартизированный коэффициент.

И res.coef выглядит следующим образом:

res.coef 
    model observ   raw   t  abst   std 
m1v1  1  V1 -0.12884211 -1.2438295 1.2438295 -0.14110975 
m2v2  2  V3 0.09258638 0.8922776 0.8922776 0.10161095 
m3v3  3  V2 0.01534989 0.1420060 0.1420060 0.01623527 

Один способ сравнить модели, чтобы построить диаграмму рассеяния всех против переменной

colnames(res.coef)[colnames(res.coef)== "variable"] <- "observ" 
dat <- melt(res.coef) 
ggplot(dat) + 
    geom_point(aes(observ,value,color=variable),size=5) + 
    theme_bw()