Análise de Correspondência - Renda x Nº de filhos - ANACOR

---

# install.packages("FactoMineR", "gglot2", "readxl")
# install.packages("devtools")
# devtools::install_github("kassambara/factoextra")
rm(list=ls(all=TRUE))
library("ggplot2")
library("FactoMineR")
library("factoextra")
library("readxl")
library("gplots")
library("corrplot")
library("graphics")
library("foreign")
library("readxl")

Tabela de contingência

# entrada de dados já tabulados 
# o comum é usar o comando table para cruzar os dados
tab <- matrix(data = c(15,27,50,43,25,37,12,8,8,13,9,10), nrow = 3, ncol = 4, byrow = T)
rownames(tab) <- c("menos de 2000", "2000 a 5000", "5000 ou mais")
colnames(tab) <- c("0", "1", "2", ">2")
tab

               0  1  2 >2
menos de 2000 15 27 50 43
2000 a 5000   25 37 12  8
5000 ou mais   8 13  9 10

# 1. convert the data as a table
dt <- as.table(as.matrix(tab))
# 2. Graph
balloonplot(t(dt), main ="Tabela de contingência - Renda x nº de Filhos", xlab ="", ylab="",
            label = T, show.margins = T)

As células com os resíduos padronizados absolutos mais elevados contribuem mais para a pontuação total do Qui-quadrado.

chisq <- chisq.test(tab)
chisq

    Pearson's Chi-squared test

data:  tab
X-squared = 41.391, df = 6, p-value = 2.425e-07

# Observed counts
chisq$observed

               0  1  2 >2
menos de 2000 15 27 50 43
2000 a 5000   25 37 12  8
5000 ou mais   8 13  9 10

# Expected counts 
round(chisq$expected,2)

                  0     1     2    >2
menos de 2000 25.21 40.45 37.30 32.04
2000 a 5000   15.32 24.57 22.65 19.46
5000 ou mais   7.47 11.98 11.05  9.49

# residuals padronizados 
round(chisq$residuals, 2)

                  0     1     2    >2
menos de 2000 -2.03 -2.11  2.08  1.94
2000 a 5000    2.47  2.51 -2.24 -2.60
5000 ou mais   0.19  0.29 -0.62  0.16

corrplot(chisq$residuals, is.cor = FALSE)

A inércia total (\(\phi^2\)) é a quantidade de informação contida na tabela de dados

phi2 <- as.numeric(chisq$statistic/sum(tab))
phi2

[1] 0.161053

Análise de correspondência

Análise de correspondência (ANACOR) é necessária para a tabela de contingência grande.

Aplica-se para visualizar pontos da linha e pontos da coluna em um espaço dimensional reduzido.

ANACOR é um método de redução dimensional aplicado a uma tabela de contingência.

A informação retida por cada dimensão é chamada autovalor.

A informação total (ou inércia) contida nos dados é chamada phi (\(\phi^2\)) e pode ser calculado da seguinte forma:

\(\phi^2 = \frac{\chi^2}{grand.total}\)

A análise de correspondência é usada para representar graficamente a tabela de distâncias entre variáveis de linha ou entre variáveis de coluna.

A abordagem ANACOR inclui as seguintes etapas:

1.Calcular os resíduos padronizados

Os resíduos padronizados (S) são:

\(S = \frac{o - e}{\sqrt{e}}\)

De fato, S é apenas a raíz quadrada dos termos que compõem a estatística \(\chi^2\).

2. Calcule a decomposição do valor singular (SVD) dos resíduos padronizados.

\(M = \frac{1}{sqrt(grand.total)} \times S\)

SVD significa que queremos encontrar matrizes ortogonais U e V, em conjunto com uma matriz diagonal \(\Delta\), tal que:

\(M = U \Delta V^T\)

• U É uma matriz contendo autovetores de linhas
• \(\Delta\) É a matriz diagonal. Os números na diagonal da matriz são chamados de valores singulares (SV). Os autovalores são o SV quadrado.
• V É uma matriz contendo autovetores de coluna

O autovalor de um determinado eixo é:

\(\lambda = \delta^2\)

\(\delta\) É o valor singular

As coordenadas das variáveis de linha em um determinado eixo são:

\(row.coord = \frac{U * \delta }{\sqrt{row.mass}}\)

As coordenadas das colunas são

\(col.coord = \frac{V * \delta }{\sqrt{col.mass}}\)

cálculo SVD

### Grand total
n <- sum(tab)
n

[1] 257

### Standardized residuals
residuals <- chisq$residuals/sqrt(n)
residuals

                        0          1           2          >2
menos de 2000 -0.12688453 -0.1318950  0.12976393  0.12074460
2000 a 5000    0.15437047  0.1564536 -0.13962533 -0.16207938
5000 ou mais   0.01207689  0.0182992 -0.03847852  0.01024036

### Number of dimensions
nb.axes <- min(nrow(residuals)-1, ncol(residuals)-1)
nb.axes

[1] 2

### Singular value decomposition
res.svd <- svd(residuals, nu = nb.axes, nv = nb.axes)
res.svd

$d
[1] 3.994798e-01 3.832625e-02 4.351618e-17

$u
            [,1]       [,2]
[1,]  0.63728969  0.2618588
[2,] -0.76722779  0.3038015
[3,] -0.07227289 -0.9160429

$v
           [,1]        [,2]
[1,] -0.5010825  0.06807878
[2,] -0.5142022 -0.09836301
[3,]  0.4821335  0.69950511
[4,]  0.5020554 -0.70454425

### singular value
sv <- res.svd$d[1:nb.axes] 
sv

[1] 0.39947979 0.03832625

u <-res.svd$u
u

            [,1]       [,2]
[1,]  0.63728969  0.2618588
[2,] -0.76722779  0.3038015
[3,] -0.07227289 -0.9160429

v <- res.svd$v
v

           [,1]        [,2]
[1,] -0.5010825  0.06807878
[2,] -0.5142022 -0.09836301
[3,]  0.4821335  0.69950511
[4,]  0.5020554 -0.70454425

# inércia total
sum(sv^2)

[1] 0.161053

### Eigenvalues
eig <- sv^2
### Variances in percentage
variance <- eig*100/sum(eig)
### Cumulative variances
cumvar <- cumsum(variance)
eig<- data.frame(eig = eig, variance = variance,
                     cumvariance = cumvar)
head(eig)

barplot(eig[, 2], names.arg=1:nrow(eig), 
       main = "Variances",
       xlab = "Dimensions",
       ylab = "Percentage of variances",
       col ="steelblue")
### Add connected line segments to the plot
lines(x = 1:nrow(eig), eig[, 2], 
      type="b", pch=19, col = "red")

Quantas dimensões reter ?

1. O número máximo de eixos na ANACOR é:

\(nb.axes = min(r-1, c-1)\)

r e c são respectivamente o número de linhas e colunas na tabela.

Coordenadas principais das linhas (Row coordinates)

### row sum
row.sum <- apply(tab, 1, sum)
row.sum

menos de 2000   2000 a 5000  5000 ou mais 
          135            82            40 

### row mass
row.mass <- row.sum/n
row.mass

menos de 2000   2000 a 5000  5000 ou mais 
    0.5252918     0.3190661     0.1556420 

### row coord = sv * u /sqrt(row.mass)
cc <- t(apply(u, 1, '*', sv)) ### each row X sv
row.coord <- apply(cc, 2, '/', sqrt(row.mass))
rownames(row.coord) <- rownames(tab)
colnames(row.coord) <- paste0("Dim.", 1:nb.axes)
round(row.coord,3)

               Dim.1  Dim.2
menos de 2000  0.351  0.014
2000 a 5000   -0.543  0.021
5000 ou mais  -0.073 -0.089

### plot
plot(row.coord, pch=19, col = "blue")
text(row.coord, labels =rownames(row.coord), pos = 3, col ="blue")

abline(v=0, h=0, lty = 2)

Coordenadas principais das colunas (Column coordinates)

### Coordinates of columns
col.sum <- apply(tab, 2, sum)
col.mass <- col.sum/n
### coordinates sv * v /sqrt(col.mass)
cc <- t(apply(v, 1, '*', sv))
col.coord <- apply(cc, 2, '/', sqrt(col.mass))
rownames(col.coord) <- colnames(tab)
colnames(col.coord) <- paste0("Dim", 1:nb.axes)
head(col.coord)

         Dim1         Dim2
0  -0.4631802  0.006037457
1  -0.3752753 -0.006887306
2   0.3664372  0.051006395
>2  0.4116690 -0.055425073

### plot
plot(col.coord, pch=17, col = "red")
text(col.coord, labels =rownames(col.coord), pos = 3, col ="red")

abline(v=0, h=0, lty = 2)

Biplot de linhas e colunas para ver a associação

xlim <- range(c(row.coord[,1], col.coord[,1]))*1.1
ylim <- range(c(row.coord[,2], col.coord[,2]))*1.1
### Plot of rows
plot(row.coord, pch=19, col = "blue", xlim = xlim, ylim = ylim)
text(row.coord, labels =rownames(row.coord), pos = 3, col ="blue")

### plot off columns
points(col.coord, pch=17, col = "red")
text(col.coord, labels =rownames(col.coord), pos = 3, col ="red")

abline(v=0, h=0, lty = 2)

Diagnóstico

Lembre-se que, a inércia total contida nos dados é:

\(\phi^2 = \frac{\chi^2}{n} = 0.3176\)

Nosso gráfico bidimensional captura 100% da inércia total da tabela.

Contribuição de linhas e colunas

As contribuições de uma linha / coluna para a definição de um eixo principal são:

\(row.contrib = \frac{row.mass * row.coord^2}{eigenvalue}\)

\(col.contrib = \frac{col.mass * col.coord^2}{eigenvalue}\)

Contribuição de linhas em %

### contrib <- row.mass * row.coord^2/eigenvalue
cc <- apply(row.coord^2, 2, "*", row.mass)
row.contrib <- t(apply(cc, 1, "/", eig[1:nb.axes,1])) *100
round(row.contrib, 2)

              Dim.1 Dim.2
menos de 2000 40.61  6.86
2000 a 5000   58.86  9.23
5000 ou mais   0.52 83.91

corrplot(row.contrib, is.cor = FALSE)

Contribuição das colunas em %

### contrib <- col.mass * col.coord^2/eigenvalue
cc <- apply(col.coord^2, 2, "*", col.mass)
col.contrib <- t(apply(cc, 1, "/", eig[1:nb.axes,1])) *100
round(col.contrib, 2)

    Dim1  Dim2
0  25.11  0.46
1  26.44  0.97
2  23.25 48.93
>2 25.21 49.64

corrplot(col.contrib, is.cor = FALSE)

Qualidade da representação

A qualidade da representação é chamada COS2.

A qualidade da representação de uma linha em um eixo é:

\(row.cos2 = \frac{row.coord^2}{d^2}\)

• Row.coord é a coordenada da linha no eixo
• \(d^2\) É a distância ao quadrado do perfil médio

Lembre-se de que a distância entre cada perfil de linha e o perfil de linha médio é:

\(d^2(row_i, average.profile) = \sum{\frac{(row.profile_i - average.profile)^2}{average.profile}}\)

row.profile <- tab/row.sum
head(round(row.profile, 3))

                  0     1     2    >2
menos de 2000 0.111 0.200 0.370 0.319
2000 a 5000   0.305 0.451 0.146 0.098
5000 ou mais  0.200 0.325 0.225 0.250

average.profile <- col.sum/n
head(round(average.profile, 3))

    0     1     2    >2 
0.187 0.300 0.276 0.237 

O código R abaixo calcula a distância do perfil médio para todas as variáveis de linha

average.rp <- col.sum/n 
d2.row <- apply(row.profile, 1, 
                function(row.p, av.p){sum(((row.p - av.p)^2)/av.p)}, 
                average.rp)
head(round(d2.row,3))

menos de 2000   2000 a 5000  5000 ou mais 
        0.124         0.295         0.013 

O cos2 de linhas no mapa de fatores são:

row.cos2 <- apply(row.coord^2, 2, "/", d2.row)
round(row.cos2, 3)

              Dim.1 Dim.2
menos de 2000 0.998 0.002
2000 a 5000   0.999 0.001
5000 ou mais  0.403 0.597

corrplot(row.cos2, is.cor = FALSE)

Cos2 das colunas

\(col.cos2 = \frac{col.coord^2}{d^2}\)

col.profile <- t(tab)/col.sum
col.profile <- t(col.profile)
###head(round(col.profile, 3))
average.profile <- row.sum/n
###head(round(average.profile, 3))

O código R abaixo calcula a distância do perfil médio para todas as variáveis da coluna

d2.col <- apply(col.profile, 2, 
        function(col.p, av.p){sum(((col.p - av.p)^2)/av.p)}, 
        average.profile)
round(d2.col,3)

    0     1     2    >2 
0.215 0.141 0.137 0.173 

O cos2 das colunas no mapa de fatores são:

col.cos2 <- apply(col.coord^2, 2, "/", d2.col)
round(col.cos2, 3)

    Dim1  Dim2
0  1.000 0.000
1  1.000 0.000
2  0.981 0.019
>2 0.982 0.018

corrplot(col.cos2, is.cor = FALSE)

Pacotes no R para ANACOR

-FactoMineR -ade4 -ca

library(FactoMineR)
res.ca <- CA(tab, graph = F)
summary(res.ca)

Call:
CA(X = tab, graph = F) 

The chi square of independence between the two variables is equal to 41.39062 (p-value =  2.425291e-07 ).

Eigenvalues
                       Dim.1   Dim.2
Variance               0.160   0.001
% of var.             99.088   0.912
Cumulative % of var.  99.088 100.000

Rows
                Iner*1000    Dim.1    ctr   cos2    Dim.2    ctr   cos2  
menos de 2000 |    64.914 | -0.351 40.614  0.998 | -0.014  6.857  0.002 |
2000 a 5000   |    94.073 |  0.543 58.864  0.999 | -0.021  9.230  0.001 |
5000 ou mais  |     2.066 |  0.073  0.522  0.403 |  0.089 83.913  0.597 |

Columns
                Iner*1000    Dim.1    ctr   cos2    Dim.2    ctr   cos2  
0             |    40.076 |  0.463 25.108  1.000 | -0.006  0.463  0.000 |
1             |    42.209 |  0.375 26.440  1.000 |  0.007  0.968  0.000 |
2             |    37.815 | -0.366 23.245  0.981 | -0.051 48.931  0.019 |
>2            |    40.954 | -0.412 25.206  0.982 |  0.055 49.638  0.018 |

### eigenvalue
head(res.ca$eig)[, 1:2]

### barplot of percentage of variance
barplot(res.ca$eig[,2], names.arg = rownames(res.ca$eig))

### Plot row points
plot(res.ca, invisible ="col")

### Plot column points
plot(res.ca, invisible ="row")

### Biplot of rows and columns
plot(res.ca)

plot(res.ca)

ellipseCA(res.ca,ellipse="col")

ellipseCA(res.ca,ellipse="row")

ellipseCA(res.ca)