Há algum tempo venho analizando um grupo de dados coletados em uma lagoa do Litoral do Rio Grande do Sul com o intuito de modelar a distribuição de moluscos limnicos. Venho usando o software R para estas análises. Um dia posto aqui os prós e os contras de usar o programa R, apesar de muita gente achar que o R é impossível.

Como uma primeira abordagem decidi investigar a colinearidade das variáveis que disponho para o ajuste de um modelo com o menor número possível de variáveis. Dando uma buscada no blog R-Bloggers descrobri uma função criada por Stephen Turner (http://gettinggeneticsdone.blogspot.com.br/2012/08/more-on-exploring-correlations-in-r.html) que coloca informações relevantes a cerca dos dados que disponho como o coeficiente de correlação e o valor de p(asteriscos) em um formato de matriz com todas as variáveis que quero analizar. Bem útil para ter uma visão da colinearidade dos preditores e o uso das funções são simples, sem problemas. Esse é um exemplo que o ambiente R e um ferramenta para a criatividade. Para mais detalhes acesse o post original. Abaixo reproduzo o código das funções para serem rodadas no R, usadas e aprovadas.

Função que cria a matriz de correlação com o valor de p.

cor.prob <- function(X, dfr = nrow(X) - 2) {
    R <- cor(X, use = "pairwise.complete.obs")
    above <- row(R) < col(R)
    r2 <- R[above]^2
    Fstat <- r2 * dfr/(1 - r2)
    R[above] <- 1 - pf(Fstat, 1, dfr)
    R[row(R) == col(R)] <- NA
    R
}

Essa função “achata” a matriz criando uma matriz com quatro colunas: uma coluna para a variável i, uma para a variável j, uma para a correlação das variáveis e uma última coluna para o valor de p.

flattenSquareMatrix <- function(m) {
    if ((class(m) != "matrix") | (nrow(m) != ncol(m))) 
        stop("Must be a square matrix.")
    if (!identical(rownames(m), colnames(m))) 
        stop("Row and column names must be equal.")
    ut <- upper.tri(m)
    data.frame(i = rownames(m)[row(m)[ut]], j = rownames(m)[col(m)[ut]], cor = t(m)[ut], 
        p = m[ut])
}

Agora selecione uma matriz de dados, aqui o autor usou uma matriz de dados mtcars do R.

mydata <- mtcars[, c(1, 3, 4, 5, 6)]

Gerando a matriz de correlação com a função cor().

cor(mydata)
##          mpg    disp      hp    drat      wt
## mpg   1.0000 -0.8476 -0.7762  0.6812 -0.8677
## disp -0.8476  1.0000  0.7909 -0.7102  0.8880
## hp   -0.7762  0.7909  1.0000 -0.4488  0.6587
## drat  0.6812 -0.7102 -0.4488  1.0000 -0.7124
## wt   -0.8677  0.8880  0.6587 -0.7124  1.0000

Matriz de correlação com o valor de p.

cor.prob(mydata)
##          mpg       disp         hp       drat        wt
## mpg       NA  9.380e-10  1.788e-07  1.776e-05 1.294e-10
## disp -0.8476         NA  7.143e-08  5.282e-06 1.222e-11
## hp   -0.7762  7.909e-01         NA  9.989e-03 4.146e-05
## drat  0.6812 -7.102e-01 -4.488e-01         NA 4.784e-06
## wt   -0.8677  8.880e-01  6.587e-01 -7.124e-01        NA

“Achatando” a matriz.

flattenSquareMatrix(cor.prob(mydata))
##       i    j     cor         p
## 1   mpg disp -0.8476 9.380e-10
## 2   mpg   hp -0.7762 1.788e-07
## 3  disp   hp  0.7909 7.143e-08
## 4   mpg drat  0.6812 1.776e-05
## 5  disp drat -0.7102 5.282e-06
## 6    hp drat -0.4488 9.989e-03
## 7   mpg   wt -0.8677 1.294e-10
## 8  disp   wt  0.8880 1.222e-11
## 9    hp   wt  0.6587 4.146e-05
## 10 drat   wt -0.7124 4.784e-06

Gerando o lindo gráfico, genial.

library(PerformanceAnalytics)
## Loading required package: zoo
## Attaching package: 'zoo'
## The following object(s) are masked from 'package:base':
## 
## as.Date, as.Date.numeric
## Loading required package: xts
## Attaching package: 'PerformanceAnalytics'
## The following object(s) are masked from 'package:graphics':
## 
## legend
chart.Correlation(mydata)

plot of chunk unnamed-chunk-7