Como reorganizar dados 2D para obter correlação?

Eu tenho o seguinte conjunto de dados simples com duas variáveis ​​contínuas; ou seja:

d = data.frame(x=runif(100,0,100),y = runif(100,0,100))
plot(d$x,d$y)
abline(lm(y~x,d), col="red")
cor(d$x,d$y) # = 0.2135273

Preciso reorganizar os dados para que a correlação entre as variáveis ​​seja ~ 0,6. Eu preciso manter os meios e outras estatísticas descritivas (sd, min, max, etc.) De ambas as variáveis ​​constantes.

Eu sei que é possível fazer quase qualquer correlação com os dados fornecidos, ou seja:

d2 = with(d,data.frame(x=sort(x),y=sort(y)))
plot(d2$x,d2$y)
abline(lm(y~x,d2), col="red")
cor(d2$x,d2$y) # i.e. 0.9965585

Se eu tentar usar a samplefunção para esta tarefa:

cor.results = c()
for(i in 1:1000){
    set.seed(i)
    d3 = with(d,data.frame(x=sample(x),y=sample(y)))
    cor.results =  c(cor.results,cor(d3$x,d3$y))
}

Eu recebo uma ampla gama de correlações:

> summary(cor.results)
     Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
-0.281600 -0.038330 -0.002498 -0.001506  0.034380  0.288800

mas esse intervalo depende do número de linhas no quadro de dados e da diminuição com o aumento do tamanho.

> d = data.frame(x=runif(1000,0,100),y = runif(1000,0,100))
> cor.results = c()
> for(i in 1:1000){
+ set.seed(i)
+ d3 = with(d,data.frame(x=sample(x),y=sample(y)))
+ cor.results =  c(cor.results,cor(d3$x,d3$y))
+ }
> summary(cor.results)
      Min.    1st Qu.     Median       Mean    3rd Qu.       Max. 
-0.1030000 -0.0231300 -0.0005248 -0.0005547  0.0207000  0.1095000

Minha pergunta é:

Como reorganizar esse conjunto de dados para obter uma correlação determinada (ou seja, 0,7)? (Também será bom se o método remover a dependência do tamanho do conjunto de dados)

Aqui está uma maneira de reorganizar os dados com base na geração de números aleatórios adicionais.

Extraímos amostras de uma distribuição normal bivariada com correlação especificada. Em seguida, calculamos as classificações dos valores e que obtemos. Essas classificações são usadas para ordenar os valores originais. Para essa abordagem, classificamos os valores e originais .$x$$y$$x$$y$

Primeiro, criamos o conjunto de dados real (como no seu exemplo).

set.seed(1)
d <- data.frame(x = runif(100, 0, 100), y = runif(100, 0, 100))

cor(d$x, d$y)
# [1] 0.01703215

Agora, especificamos uma matriz de correlação.

corr <- 0.7  # target correlation
corr_mat <- matrix(corr, ncol = 2, nrow = 2)
diag(corr_mat) <- 1
corr_mat
#      [,1] [,2]
# [1,]  1.0  0.7
# [2,]  0.7  1.0

Geramos dados aleatórios após uma distribuição normal bivariada com , (para ambas as variáveis) e a correlação especificada. Em R, isso pode ser feito com a função do pacote. Usamos para indicar que a correlação é a correlação empírica (não a correlação da população).$\mu = 0$$\sigma = 1$mvrnormMASSempirical = TRUE

library(MASS)
mvdat <- mvrnorm(n = nrow(d), mu = c(0, 0), Sigma = corr_mat, empirical = TRUE)

cor(mvdat)
#      [,1] [,2]
# [1,]  1.0  0.7
# [2,]  0.7  1.0

Os dados aleatórios correspondem perfeitamente à correlação especificada.

Em seguida, calculamos as classificações dos dados aleatórios.

rx <- rank(mvdat[ , 1], ties.method = "first")
ry <- rank(mvdat[ , 2], ties.method = "first")

Para usar as classificações para os dados originais d, precisamos classificar os dados originais.

dx_sorted <- sort(d$x)
dy_sorted <- sort(d$y)

Agora, podemos usar as classificações para especificar a ordem dos dados classificados.

cor(dx_sorted[rx], dy_sorted[ry])
# [1] 0.6868986

A correlação obtida não corresponde perfeitamente à especificada, mas a diferença é relativamente pequena.

Aqui, dx_sorted[rx]e dy_sorted[ry]são versões reamostradas dos dados originais em d.

Para gerar duas distribuições uniformes com uma correlação especificada, o algoritmo de Ruscio & Kaczetow (2008) funcionará. Eles fornecem código R . Você pode transformar com uma função linear simples para obter seu alvo mínimo, máximo, média e DP.

Algoritmo de Ruscio e Kaczetow

Resumirei o caso bivariado, mas também pode funcionar com problemas multivariados. Uncorrelated e são gerados com qualquer forma (por exemplo, uniforme). Em seguida, e são gerados como normais bivariados com uma correlação intermediária. e são substituídos por e maneira a preservar a classificação. Ajuste a correlação intermediária para ser maior ou menor, dependendo de r ( ) estar muito baixo ou muito alto. e são gerados como normais bivariados com a nova correlação intermediária. Repetir.$X_o$$Y_o$$X_1$$Y_1$$X_1$$Y_1$$X_0$$Y_0$$X_1,Y_1$$X_2$$Y_2$

Observe que isso é muito semelhante à solução de @Sven Hohenstein, exceto que é iterativa, portanto a correlação intermediária se aproxima cada vez mais da correlação de destino até que sejam indistinguíveis. Além disso, observe que esse algoritmo pode ser usado para gerar uma grande população (por exemplo, N = 1 milhão) a partir da qual é possível coletar amostras menores - isso é útil se você precisar de um erro de amostragem.

Para um post relacionado: Correlação e distribuições não normais

Preservando a estatística descritiva

Não há garantia de que o algoritmo produza exatamente os mesmos descritivos. No entanto, como a média e o DP de uma distribuição uniforme são determinados por seu mínimo e máximo, você pode simplesmente ajustar o mínimo e o máximo para corrigir tudo.

Deixe e serem suas variáveis ​​geradas a partir da última iteração do algoritmo Ruscio & Kaczetow, e sejam suas variáveis ​​finais que você espera ter (com descritivos de destino) e e sejam suas variáveis ​​originais em seu conjunto de dados. $X_g$$Y_g$$X_f$$Y_f$$X$$Y$

Calcule $X_f=(X_g - min(X))*(max(X)-min(x))/(max(X_g)-min(X_g))$

Faça o mesmo para$Y_f$

Referência:

Ruscio, J. & Kaczetow, W. (2008). Simulação de dados não normais multivariados usando um algoritmo iterativo. Pesquisa Comportamental Multivariada, 43, 355-381. doi: 10.1080 / 00273170802285693

Estou supondo que quando você diz "reamostrar" você quer dizer "simular", o que é mais geral. A seguir, é a maneira mais concisa que conheço para simular dados normais e bivariados com uma correlação especificada. Substitua seus próprios valores desejados por r e n.

r = .6
n = 1000
x = rnorm(n) 
z = rnorm(n) 
y = (r/(1-r^2)^.5)*x + z

cor(x,y)
plot(x,y)
abline(lm(y~x), col="red")