• Para realizar análise de dados é preciso ter os dados, mais ainda, é preciso que o dado obtido seja uma medida da sua variável de interesse. Nem sempre encontramos esses dados, fica pior, nem sempre é possível medir o que queremos.
• Aplicar enquetes é uma forma popular de obter dados para uma pesquisa, mas para isso as perguntas devem ser adequadas para capturar a variável que desejamos estudar e é preciso considerar questões de comportamento de quem responde os questionários que podem induzir um viés na análise.
• Nesta unidade vamos discutir questões relacionadas a dados obtidos por enquetes (dados de survey). Como fizemos nas outras unidades, começamos com um exemplo.

• Na sequência dos atentados de 11 de Setembro de 2001 os Estados Unidos invadiram o Afeganistão com o objetivo de destruir a al-Qaeda. Em 2003 a OTAN entrou na guerra e mandou tropas com o nome de Força Internacional de Apoio à Segurança (ISAF, do inglês International Security Assistance Force).
• A ISAF adotou uma campanha para conquistar corações e mentes da população, a campanha combinava assistência econômica, prestação de serviços e proteção aos civis.
• Para avaliar a campanha é essencial entender e medir as experiências e sentimentos dos civis durante o período de guerra. Não bastasse as dificuldades inatas de medir sentimentos, as condições de guerra podiam tornar perigoso revelar os sentimentos e as experiencias respondendo as questões das pesquisas.

• Um grupo de cientistas sociais resolveu fazer uma pesquisa de opinião no sul do Afeganistão, o coração da guerra de insurgência.
• São dois os artigos que fundamentam esse exemplo:
  • Jason Lyall, Graemer Blair e Kosuke Imai. Explaining support for combatants during wartime: a survey experiment in Afghanistan. American Political Science Review, v.107, 2013.
  • Graeme Blair, Kosuke Imai e Jason Lyall. Comparing and combining list and endorsement experiments: evidence from Afghanistan. American Journal of Political Science, v.58, 2014.

• A pesquisa foi aplicada em uma amostra de 2754 cidadãos do Afeganistão entre janeiro e fevereiro de 2011.
• A taxa de participação foi de 89%, ou seja, foram contactados 3097 homens, mas 343 se recusaram a responder.
• Por questões locais, apenas homens responderam às pesquisas.

• Variáveis:
  • province: província onde vive quem respondeu.
  • district: distrito onde vive quem respondeu.
  • village.id: o ID da vila onde vive quem respondeu.
  • age: idade.
  • educ.years: anos de educação.

• Variáveis:
  • employed: se quem respondeu está empregado.
  • income: renda mensal (cinco níveis)
  • violent.ex.ISAF: se experimentou violência pelo ISAF
  • violent.ex.taliban: se experimentou violência pelo taliban
  • list.group: grupo designado aleatoriamente para list experiment (control, ISAF, taliban)
  • list.response: repostas para as questões do list experiment (0-4)

afghan <- read.csv("afghan.csv")
str(afghan)

## 'data.frame':    2754 obs. of  11 variables:
##  $ province           : chr  "Logar" "Logar" "Logar" "Logar" ...
##  $ district           : chr  "Baraki Barak" "Baraki Barak" "Baraki Barak" "Baraki Barak" ...
##  $ village.id         : int  80 80 80 80 80 80 80 80 80 35 ...
##  $ age                : int  26 49 60 34 21 18 42 39 20 18 ...
##  $ educ.years         : int  10 3 0 14 12 10 6 12 5 10 ...
##  $ employed           : int  0 1 1 1 1 1 1 1 1 0 ...
##  $ income             : chr  "2,001-10,000" "2,001-10,000" "2,001-10,000" "2,001-10,000" ...
##  $ violent.exp.ISAF   : int  0 0 1 0 0 0 0 0 0 0 ...
##  $ violent.exp.taliban: int  0 0 0 0 0 0 0 1 0 0 ...
##  $ list.group         : chr  "control" "control" "control" "ISAF" ...
##  $ list.response      : int  0 1 1 3 3 2 1 3 1 2 ...

head(afghan)

##   province     district village.id age educ.years employed       income
## 1    Logar Baraki Barak         80  26         10        0 2,001-10,000
## 2    Logar Baraki Barak         80  49          3        1 2,001-10,000
## 3    Logar Baraki Barak         80  60          0        1 2,001-10,000
## 4    Logar Baraki Barak         80  34         14        1 2,001-10,000
## 5    Logar Baraki Barak         80  21         12        1 2,001-10,000
## 6    Logar Baraki Barak         80  18         10        1         <NA>
##   violent.exp.ISAF violent.exp.taliban list.group list.response
## 1                0                   0    control             0
## 2                0                   0    control             1
## 3                1                   0    control             1
## 4                0                   0       ISAF             3
## 5                0                   0       ISAF             3
## 6                0                   0       ISAF             2

library(tidyverse)
afghan %>%
  select(age, educ.years, employed, income) %>%
  summary()

##       age          educ.years        employed         income         
##  Min.   :15.00   Min.   : 0.000   Min.   :0.0000   Length:2754       
##  1st Qu.:22.00   1st Qu.: 0.000   1st Qu.:0.0000   Class :character  
##  Median :30.00   Median : 1.000   Median :1.0000   Mode  :character  
##  Mean   :32.39   Mean   : 4.002   Mean   :0.5828                     
##  3rd Qu.:40.00   3rd Qu.: 8.000   3rd Qu.:1.0000                     
##  Max.   :80.00   Max.   :18.000   Max.   :1.0000

• Eu definiria income como fator…

afghan %>%
  select(age, educ.years, employed, income) %>%
  mutate(income = as.factor(income)) %>%
  summary()

##       age          educ.years        employed                  income    
##  Min.   :15.00   Min.   : 0.000   Min.   :0.0000   10,001-20,000  : 616  
##  1st Qu.:22.00   1st Qu.: 0.000   1st Qu.:0.0000   2,001-10,000   :1420  
##  Median :30.00   Median : 1.000   Median :1.0000   20,001-30,000  :  93  
##  Mean   :32.39   Mean   : 4.002   Mean   :0.5828   less than 2,000: 457  
##  3rd Qu.:40.00   3rd Qu.: 8.000   3rd Qu.:1.0000   over 30,000    :  14  
##  Max.   :80.00   Max.   :18.000   Max.   :1.0000   NA's           : 154

• … mas vou seguir o livro.

count(afghan, income)

##            income    n
## 1   10,001-20,000  616
## 2    2,001-10,000 1420
## 3   20,001-30,000   93
## 4 less than 2,000  457
## 5     over 30,000   14
## 6            <NA>  154

unique(afghan$income)

## [1] "2,001-10,000"    NA                "10,001-20,000"   "less than 2,000"
## [5] "20,001-30,000"   "over 30,000"

• A média de idade dos que responderam as questões é de 32 anos.
• A grande maioria tem pouca educação (a mediana é 1)
• Cerca de 60% estavam desempregados.
• A maioria dos que responderam tem renda mensal inferior a 10.000 afghani, que na época era equivalente a 200 dólares.

• Para avaliar como a violência atingiu a população civil, foi feita a seguinte pergunta: No último ano, você ou alguém de sua família foi atingido por conta da ação de forças entrangeiras ou do taliban?
• Pergunta orginal: Over the past year, have you or anyone in your family suffered harm due to the actions of the Foreign Forces/the Taliban?

• Na pesquisa foi registrado que harm, que traduzi como atingido, se refere a injúrias físicas bem como danos ao patrimônio.
• Vamos calcular a proporção dos que foram atingidos por ações da ISAF e/ou do Taliban.

harm_props <- afghan %>%
  group_by(violent.exp.ISAF, violent.exp.taliban) %>%
  count() %>%
  ungroup() %>%
  mutate(prop = n/sum(n))

## # A tibble: 9 × 4
##   violent.exp.ISAF violent.exp.taliban     n    prop
##              <int>               <int> <int>   <dbl>
## 1                0                   0  1330 0.483  
## 2                0                   1   354 0.129  
## 3                0                  NA    22 0.00799
## 4                1                   0   475 0.172  
## 5                1                   1   526 0.191  
## 6                1                  NA    22 0.00799
## 7               NA                   0     7 0.00254
## 8               NA                   1     8 0.00290
## 9               NA                  NA    10 0.00363

• Sem usar ungroup() a soma da última linha seria feita dentro de cada grupo.

afghan %>%
  group_by(violent.exp.ISAF, violent.exp.taliban) %>%
  count() %>%
#  ungroup() %>%
  mutate(prop = n/sum(n))

## # A tibble: 9 × 4
## # Groups:   violent.exp.ISAF, violent.exp.taliban [9]
##   violent.exp.ISAF violent.exp.taliban     n  prop
##              <int>               <int> <int> <dbl>
## 1                0                   0  1330     1
## 2                0                   1   354     1
## 3                0                  NA    22     1
## 4                1                   0   475     1
## 5                1                   1   526     1
## 6                1                  NA    22     1
## 7               NA                   0     7     1
## 8               NA                   1     8     1
## 9               NA                  NA    10     1

• Proporção atingida pela ISAF.

ISAF_harm_prop <- harm_props %>%
  filter(violent.exp.ISAF == 1) %>%
  summarize(harm_prop = sum(prop)) %>%
  pull()

ISAF_harm_prop

## [1] 0.3714597

• Proporção atingida pelo taliban

talib_harm_prop <- harm_props %>%
  filter(violent.exp.taliban == 1) %>%
  summarize(harm_prop = sum(prop)) %>%
  pull()

talib_harm_prop

## [1] 0.3224401

• Proporção atingida pelo taliban e pela ISAF

both_harm_prop <- harm_props %>%
  filter(violent.exp.taliban == 1 & violent.exp.ISAF == 1) %>%
  summarize(harm_prop = sum(prop)) %>%
  pull()

both_harm_prop

## [1] 0.1909949

• Em muitos questionários aparecem perguntas não respondidas, mesmo em bases com dados observados é comum que alguns dados estejam faltando.
• Quando o dado está faltando o R registra como NA.

afghan %>%
  select(income) %>%
  slice(1:10)

##           income
## 1   2,001-10,000
## 2   2,001-10,000
## 3   2,001-10,000
## 4   2,001-10,000
## 5   2,001-10,000
## 6           <NA>
## 7  10,001-20,000
## 8   2,001-10,000
## 9   2,001-10,000
## 10          <NA>

• É fundamental ter certeza que todos os dados ausentes são reconhecidos pelo R como ausentes e registrados como NA. Algumas vezes garantir que isso aconteça é trabalhoso, algumas bases registram NAs como espaços em branco, -, –, ou números estranhos como -999. Um NA não reconhecido pelo R pode comprometer toda a análise de dados.
• A função is.na() identifica os NAs presentes nos dados.

afghan %>%
  select(income) %>%
  slice(1:10) %>%
  is.na()

##       income
##  [1,]  FALSE
##  [2,]  FALSE
##  [3,]  FALSE
##  [4,]  FALSE
##  [5,]  FALSE
##  [6,]   TRUE
##  [7,]  FALSE
##  [8,]  FALSE
##  [9,]  FALSE
## [10,]   TRUE

• Como FALSE equivale a 0 e TRUE equivale a 1, a soma dos resultados da função is.na() é o número de NAs e a média é a proporção de NAs.

afghan %>%
  summarize(n_missing = sum(is.na(income)),
            p_missing = mean(is.na(income)))

##   n_missing  p_missing
## 1       154 0.05591866

• No R os NAs são contagiosos, ou seja, operações com NAs retornam NAs.

5 + NA

## [1] NA

5*NA

## [1] NA

• Por serem contagiosos o pesquisador é forçado a decidir o que fazer com os NAs. Em alguns casos ignorá-los é seguro, em outros não. Cabe ao pesquisador decidir o que fazer.
• Por exemplo, suponha dados de estações que medem a velocidade do vento e que os aparelhos não conseguem medir velocidades muito altas. Ignorar os NAs vai levar a subestimar a velocidade do vento na região.

• Muita funções do R possuem o argumento na.rm que permite ignorar NAs. Essa não é uma prática errada, mas o pesquisador deve estar ciente dos riscos ao escolher ignorar os NAs.

x <- c(1,2,3,4,NA)
sum(x)

## [1] NA

sum(x, na.rm = TRUE)

## [1] 10

x <- c(1,2,3,4,NA)
mean(x)

## [1] NA

mean(x, na.rm = TRUE)

## [1] 2.5

• A função filter() combinada com a função is.na() oeferece uma maneira rápida de tirar NAs de uma análise.

afghan %>%
  filter(!is.na(violent.exp.ISAF), !is.na(violent.exp.taliban)) %>%
  group_by(violent.exp.ISAF, violent.exp.taliban) %>%
  count() %>%
  ungroup() %>%
  mutate(prop = n/sum(n)) %>%
  arrange(prop)

## # A tibble: 4 × 4
##   violent.exp.ISAF violent.exp.taliban     n  prop
##              <int>               <int> <int> <dbl>
## 1                0                   1   354 0.132
## 2                1                   0   475 0.177
## 3                1                   1   526 0.196
## 4                0                   0  1330 0.495

• A função arrange() organizou os dados da menor para maior proporção, podemos organizar em ordem decrescente usando a função desc().

afghan %>%
  filter(!is.na(violent.exp.ISAF), !is.na(violent.exp.taliban)) %>%
  group_by(violent.exp.ISAF, violent.exp.taliban) %>%
  count() %>%
  ungroup() %>%
  mutate(prop = n/sum(n)) %>%
  arrange(desc(prop))

## # A tibble: 4 × 4
##   violent.exp.ISAF violent.exp.taliban     n  prop
##              <int>               <int> <int> <dbl>
## 1                0                   0  1330 0.495
## 2                1                   1   526 0.196
## 3                1                   0   475 0.177
## 4                0                   1   354 0.132

• As funções filter() e is.na() também poder ser usadas para identificar NAs em mais de uma variável, por exemplo, podemos encontrar quantos NAs para quem sofreu violência com ISAF ou taliban.

missing_prop <- harm_props %>%
  filter(is.na(violent.exp.ISAF) | is.na(violent.exp.taliban)) %>%
  summarize(missing_prop = sum(prop)) %>%
  pull

missing_prop

## [1] 0.02505447

• As funções na.omit(), do R básico, e drop_na(), do tidyverse excluem as observações onde aparece NA em pelo menos uma variável. Repare que as funções excluem toda a linha mesmo que esteja faltando dados para apenas uma variável.

nrow(afghan)

## [1] 2754

afghan.sub <- na.omit(afghan)
nrow(afghan.sub)

## [1] 2554

nrow(afghan)

## [1] 2754

afghan.sub2 <- drop_na(afghan)
nrow(afghan.sub2)

## [1] 2554

afghan %>% drop_na(income) %>% nrow()

## [1] 2600

• Além de estatísitcas descritivas podemos passar informações sobre a distribuição de uma variável de interesse usando gráficos.
• Na introdução vimos como usar o ggplot2 para criar gráficos por etapas, agora vamos usar esse conhecimento para passar informações a respeito de como os civis afegãos foram atingidos pela ISAF e pelo taliban.

• Para começar vamos ilustrar as respostas sobre para a pergunta sobre ter sido atingido pela ISAF, queremos comparar as repostas afirmativas, negativas e os que não responderam.
• Gráficos de barras fazem bem esse trabalho (muito cuidado com gráficos de pizza!). Na introdução usamos a função geom_col() para fazer gráficos de barras, aqui vamos usar a função geom_bar(). A função geom_col() é uma versão simplificada da geom_bar().

• Vítimas da ISAF

ggplot(data = afghan, aes(x=as.factor(violent.exp.ISAF))) +
  geom_bar(aes(y= after_stat(prop))) +
  scale_x_discrete(labels = c("No harm", "Harm", "Nonresponse")) +
  ylab("Proportion of respondents") +
  xlab("Response category") +
  ggtitle("Civilian Victimization by the ISAF")

• No livro a função aes() no arumento da função geom_bar() é chamada da seguinte forma: geom_bar(aes(y=..prop.., group=1)) o argumento ..prop.. foi substituido por after_stat(prop) para ficar compatível com as versões mais modernas do ggplot2. O argumento group = 1 é um “truque” para fazer o agrupamento pelos níveis de violent.exp.ISAF. Por fim, note que os autores transformaram a variável violent.exp.ISAF em fator.
• As barras podem ser feitas por vários critérios, o argumento prop na função after_stat() diz para usar proporções. Para ilustrar a função after_stat() vou refazer o gráfico usando o argumento count.

ggplot(data = afghan, aes(x=as.factor(violent.exp.ISAF))) +
  geom_bar(aes(y= after_stat(count))) +
  scale_x_discrete(labels = c("No harm", "Harm", "Nonresponse")) +
  ylab("Proportion of respondents") +
  xlab("Response category") +
  ggtitle("Civilian Victimization by the ISAF")

• Vítimas do taliban

ggplot(data = afghan, aes(x=as.factor(violent.exp.taliban))) +
  geom_bar(aes(y= after_stat(prop))) +
  scale_x_discrete(labels = c("No harm", "Harm", "Nonresponse")) +
  ylab("Proportion of respondents") +
  xlab("Response category") +
  ggtitle("Civilian Victimization by the ISAF")

• Os gráficos ilustram bem cada uma das distribuições, por exemplo, é fácil ver que a proporção de No harm é maior do que de haarm nas duas distribuições.
• Porém, não é fácil comparar um distribuição com a outra por meio dos gráficos. Por exemplo, fica difícil saber se houve mais harm para ISAF ou taliban apenas olhando os gráficos.
• Uma nova arrumação dos dados permite criar um gráfico onde a comparação é mais fácil.

afghan %>%
  pivot_longer(violent.exp.ISAF:violent.exp.taliban, 
               names_to ="harming group", values_to="harm") %>%
  ggplot(aes(x=as.factor(harm))) +
  geom_bar(aes(y=after_stat(prop), fill=harming_group, group=harming_group),
           position = "dodge") +
  scale_x_discrete(labels=c("No harm", "Harm", "Nonresponse")) +
  scale_fill_discrete(name="Harming group", labels=c("ISAF", "Taliban")) +
  labs(title = "Civilian victimization",
       x = "Response category",
       y = "Proportion of respondents")

• Para entender melhor o que aconteceu é válido observar o data.frame original e como foi modificado nas linhas que antecedem a função ggplot()

afghan %>%
  pivot_longer(violent.exp.ISAF:violent.exp.taliban, 
               names_to ="harming group", values_to="harm")

• Que tal mudar a paleta de cores e o tema?

afghan %>%
  pivot_longer(violent.exp.ISAF:violent.exp.taliban, 
               names_to ="harming_group", values_to="harm") %>%
  ggplot(aes(x=as.factor(harm))) +
  geom_bar(aes(y=after_stat(prop), fill=harming_group, group=harming_group),
           position = "dodge") +
  scale_x_discrete(labels=c("No harm", "Harm", "Nonresponse")) +
  scale_fill_brewer(name="Harming group", labels=c("ISAF", "Taliban"), 
                      palette = "Blues") +
  labs(title = "Civilian victimization",
       x = "Response category",
       y = "Proportion of respondents") +
  theme_classic()

• O histograma é um método popular para ilustrar distribuições de variáveis numéricas. Suponha, por exemplo, que queiramos dar uma olhada na idade de quem respondeu às enquetes. Como temos muitos valores um gráfico de barras ficaria confuso e não nos daria muita informação.

• Uma alternativa é dividir as idades em intervalos iguais (bins), então contamos quantas observações estão em cada intervalo (janela) e calculamos a densidade de acordo com a fórmula: \[\mbox{densidade}=\frac{\mbox{proporção de observações na janela}}{\mbox{tamanho da janela}}\]
• Em geral estamos menos preocupados com o valor exato da densidade do que em comparar a densidade das várias janelas.
• A função geom_histogram() do ggplot2 faz o histograma.

afghan %>%
  ggplot(aes(x = age)) +
  geom_histogram(aes(y= after_stat(density)), binwidth=5, boundary=0) +
  scale_x_continuous(breaks = seq(20,80,by=10)) +
  labs(title = "Distribution of respondent's age",
       y = "Density",
       x = "Age") +
  theme_classic()

• Em um histograma a área do retângulo associado a cada janela representa a proporção de observações na janela.
• Desta forma, em geral, interpretamos a escala da densidade, a unidade no eixo vertical, como percentagem pela unidade do eixo horizontal. No exemplo, a densidade é medida como percentual por idade.
• Por consequência, as alturas dos retângulos associados às janelas pode somar mais do que um, já as áreas dos retângulos somam um.

• Um histograma divide os dados em janelas onde a área de cada janela representa a proporção de observações que pertencem à janela. A altura de cada janela é chamada de densidade, que é igual a proporção de observações em cada janela dividida pela largura da janela. O histograma aproxima a distribuição de uma variável.

• Histograma como número de observações no lugar da densidade.

afghan %>%
  ggplot(aes(x = age)) +
  geom_histogram(aes(y= after_stat(count)), binwidth=5, boundary=0) +
  scale_x_continuous(breaks = seq(20,80,by=10)) +
  labs(title = "Distribuição das idades",
       y = "Número de observações",
       x = "Idade") +
  theme_classic()

• A função geom_vline() coloca linhas verticais em um gráfico (para linhas horizontais use a função geom_hline()) e a função annotate() faz anotações no gráfico.
• Vamos fazer o histograma dos anos de educação destacando a mediana.

afghan %>%
  ggplot(aes(x = educ.years, y= after_stat(density))) +
  geom_histogram(binwidth=1, center=0) +
  geom_vline(xintercept = median(afghan$educ.years)) +
  annotate(geom="text", x = median(afghan$educ.years),
           y=0.4, label="Median", hjust=-0.1)
  labs(title = "Distribution of respondent's education",
       y = "denssity",
       x = "Years of education") +
  theme_classic()

• Outro tipo de figura muito usado para ilustrar distribuições é o boxplot, diagrama de caixas. O boxplot é particularmente útil para comparar distribuições de diferentes variáveis.
• O boxplot desenha um retângulo (caixa) com 50% das observações. A altura da caixa é o intervalo interquartílico, a mediana é marcada por uma linha na caixa, abaixo da caixa tem um linha que mostra 1,5 IQR abaixo do primeiro quartil e acima da caixa uma linha que mostra 1,5 acima do terceiro quartil. Os pontos abaixo e acima da linha, pontos extremos ou outliers, são destacado na figura.
• Na ausência de pontos extremos baixos, a linha inferior acaba na menor observação. Da mesma forma, na ausência de pontos extremos altos, a linha superior acaba na maior observação.

• Boxplot da idade

afghan %>%
  ggplot(aes(y=age)) +
  geom_boxplot() +
  labs(y="Age", x=NULL, title = "Distribution of age")

• Como foi dito, o boxplot é particularmente útil para comparar distribuições. Por exemplo, suponha que você quer comparar a distribuição dos anos de educação nas diversas províncias onde a pesquisa foi aplicada.

afghan %>%
  ggplot(aes(y=educ.years, x=province)) +
  geom_boxplot() +
  labs(y="Years of education", x = "Province", title = "Education by province")

• É fácil ver que os participantes das províncias de Helmand e Uruzgam possuem menos anos de educação do que os das outras províncias, talvez não por coincidência sejam os que mais reportam ter sido alvo da ISAF e do taliban.

afghan %>%
  group_by(province) %>%
  summarize(violence.exp.taliban = mean(violent.exp.taliban, na.rm=TRUE),
            violence.exp.ISAF = mean(violent.exp.ISAF, na.rm=TRUE))

## # A tibble: 5 × 3
##   province violence.exp.taliban violence.exp.ISAF
##   <chr>                   <dbl>             <dbl>
## 1 Helmand                0.504              0.541
## 2 Khost                  0.233              0.242
## 3 Kunar                  0.303              0.399
## 4 Logar                  0.0802             0.144
## 5 Uruzgan                0.455              0.496

• Um diagrama de caixa ou boxplotmostra a distribuição de uma variável indicando a mediana, primeiro e terceiro quartis, e os pontos extremos. Permite uma comparação entre as distribuições de um conjunto de variáveis.

• A maneira mais fácil de salvar figuras é usando os menus do RSTudio, mas, para quem não gosta de seguir menus de forma alguma, existe a função ggsave(). Deixo por conta de vocês olhar o help desta função.

• Algumas vezes, por questões de estética ou para facilitar comprações, desejamos colocar figuras uma ao lado da outra. A função grid.arrange() do pacote gridExtra facilita esse trabalho.

library(gridExtra)

age_hist <- afghan %>%
  ggplot(aes(x=age)) +
  geom_histogram(aes(y= after_stat(density)), binwidth=5, boundary=0) +
  scale_x_continuous(breaks = seq(20,80,by=10)) +
  labs(title = "Distribution of respondent's age",
       y = "Density",
       x = "Age") +
  theme_classic()

educ_hist <- afghan %>%
  ggplot(aes(x = educ.years, y= after_stat(density))) +
  geom_histogram(binwidth=1, center=0) +
  geom_vline(xintercept = median(afghan$educ.years)) +
  annotate(geom="text", x = median(afghan$educ.years),
           y=0.4, label="Median", hjust=-0.1) +
  labs(title = "Distribution of respondent's education",
       y = "denssity",
       x = "Years of education") +
  theme_classic()

grid.arrange(age_hist, educ_hist, ncol=2)

• Pesquisa por amostragem é uma das principais ferramentas para cientistas sociais que fazem trabalhos aplicados. Normalmente são usadas para avaliar comportamento ou opinião da população quando as informações não estão disponíveis em outras fontes como, por exemplo, registros administrativos.
• Na pesquisa por amostragem o pesquisador seleciona um subconjunto da população, chamado de amostra, para entender características da população como um todo.
• Pesquisa por amostragem é diferente de censo, esse último tem como objetivo abranger toda a população.

• O uso de amostras permite obter informações sobre uma grande população entrevistando uma pequena fração dessa população.
• No exemplo do Afeganistão, uma amostra de 2754 pessoas foi usada para inferir experiências e comportamentos de uma população de mais de 15 milhões de pessoas.

• Nos EUA, amostras com pouco mais de mil pessoas costumam ser usadas para inferir a opinião de mais de 200 milhões de adultos.
• A “mágica” da amostra deve muito ao uso da escolha aleatória, aqui vamos tratar da amostragem aleatória simples.

• Antes do uso de amostras aleatórias, pesquisadores costumavam construir amostras por cotas. Grosso modo, o método consiste em determinar algumas caracteríticas da amostra para reproduzir a da população. Por exemplo, se 20% da população tem nível superior, então 20% da amostra deveria ser composta por pessoas com nível superior. Para definição das cotas é comum utilizar de multiplas variáveis.

• O problema das cotas é semelhante ao enfrentado por estudos com dados observados. Mesmo que a amostra reproduza caracterísitcas observadas da população que forma usadas para definir as cotas, é possível que a amostra tenha diferenças relevantes da população em termos de características não observadas.
• O uso de amostras aleatórias elimita esse viés de seleção da amostra garantindo que a amostra selecionada seja representativa da população alvo.

• Amostragem aleatória simples (SRS, do inglês: simple random sample) é o método mais simples de amostragem probabilística que evita problemas de viés de seleção da amostra por meio de escolha aleatória de unidades da população. No SRS, o número predeterminado de unidades é selecionado aleatoriamente e sem reposição de uma população-alvo, cada unidade tem a mesma probabilidade de ser selecionada. A amostra resultante é representativa da população em termos de características observadas e não observadas.

• O método de amostragem aleatória simples muitas vezes não é viável por conta de questões logísticas.
• Imagine que você vai aplicar esse método para fazer uma pesquisa aqui no DF. Onde conseguir uma lista com todos os moradores? Você pode ter pensado na lista telefônica, mas e quem não tem telefone? E quem tem vários telefones?
• Outras vezes pode ser muito difícil ou muito caro aplicar o questionário em alguns dos indivíduos escolhidos aleatóriamente.
• Por conta disso existe métodos mais sofisitcados de amostragem aleatórias.

• Para pesquisa do Afeganistão foi usado método de amostragem por cluster multiestágios. Esse método consiste em primeiro selecionar unidades grandes e depois escolher aleatoriamente unidades menores dentro das unidades maiores.
• No caso de Afeganistão, para cada um das cinco províncias de interesse foram selecionadas amostras de distritos e depois de vilas. Em cada vila os entrevistados foram escolhidos de maneira aproximadamente aleatória (a escolha levou em conta a localização na vila).

• O arquivo afghan-village.csv contém dados de altitude e população das vilas do Afeganistão.
• Vamos comparar essas propriedade nas vilas que entraram e não entraram na amostra, mas antes, vamos carregar e dar uma olhada nos dados.

afghan.village <- read.csv("afghan-village.csv")
head(afghan.village)

##   altitude population village.surveyed
## 1  1959.08        197                1
## 2  2425.88        744                0
## 3  2236.60        179                1
## 4  1691.76        225                0
## 5  1928.04        379                0
## 6  1194.56        617                0

• Histograma da altitude

afghan.village %>%
  ggplot(aes(x=altitude)) +
  geom_histogram(aes(y=after_stat(density))) +
  labs(x = "Altitude (metros)",
       y= "densidade")

• Histograma da população

afghan.village %>%
  ggplot(aes(x=population)) +
  geom_histogram(aes(y=after_stat(density))) +
  labs(x = "População",
       y= "densidade")

• Histograma do log da população

afghan.village %>%
  ggplot(aes(x=log(population))) +
  geom_histogram(aes(y=after_stat(density))) +
  labs(x = "População (log)",
       y= "densidade")

• Comparação entre as altitudes das vilas na amostra e fora da amostra

afghan.village %>%
  ggplot(aes(x = as.factor(x=village.surveyed),
             y = altitude)) +
  geom_boxplot() +
  scale_x_discrete(labels = c("Fora da amostra", "Na amostra")) +
  labs(y = "Altitude (metros)", x=NULL) +
  theme_classic()

• Comparação entre as populações das vilas na amostra e fora da amostra

afghan.village %>%
  ggplot(aes(x = as.factor(x=village.surveyed),
             y = log(population))) +
  geom_boxplot() +
  scale_x_discrete(labels = c("Fora da amostra", "Na amostra")) +
  labs(y = "População (log)", x=NULL) +
  theme_classic()

• Muitas vezes o cientista social é desafiado a desenvolver medidas para descrever ou entender comportamentos, atitudes e caracterísitcas não observadas das pessoas. Por exemplo, quantificar a ideologia de atores políticos como legisladores e juízes a partir do comportamento observado desses atores.
• Como medir o quão liberal ou conservador (lembre que o livro texto é “americano”) é um político ou, talvez mais próximo da realidade brasileira, quão de esquerda ou de direita é um senador?

• Uma estratégia adotada por cientistas políticos é observar os votos do congressista.
• O modelo de “espaço de votações”, spatial voting model, relacicona ideologia e votos de congressistas. O modelo apresenta duas dimensões: liberalismo/conservadorismo econômico e racial (pesquisadores identificaram essas como as princípais características ideológicas nos EUA do pós-guerra).
• Pesquisas mostram que o liberalismo/conservadorismo econômico tem mais poder explicativo para o voto do que o racial.

• De acordo com o modelo o congressista tem um ponto ideal no plano definido pelas dimensões relevantes, no caso liberalismo/conservadorismo racial e econômico.
• O congressista está mais propenso a votar contra uma proposta se o ponto ideal dele for mais próximo do status quo do que da proposta, vale o contrário.
• Propostas controversas revelam muito sobre a ideologia do congressista e propostas unânimes não nos dão informações sobre a ideologia do congressista.

• O exemplo de “medir ideologia” será usado para ilustrar relações bivariadas.
• A referência para o exemplo é: Nolan McCarty, Keiht T. Poole e Howard Rosenthal (2006). Polarized America: The Dance of Ideology and Unequal Richies. MIT Press.

congress <- read.csv("congress.csv")
head(congress)

##   congress district   state    party        name dwnom1 dwnom2
## 1       80        0     USA Democrat      TRUMAN -0.276  0.016
## 2       80        1 ALABAMA Democrat  BOYKIN  F. -0.026  0.796
## 3       80        2 ALABAMA Democrat   GRANT  G. -0.042  0.999
## 4       80        3 ALABAMA Democrat ANDREWS  G. -0.008  1.005
## 5       80        4 ALABAMA Democrat   HOBBS  S. -0.082  1.066
## 6       80        5 ALABAMA Democrat   RAINS  A. -0.170  0.870

• Variáveis:
  • name: nome do congressista.
  • state: estado do congressista.
  • party: partido do congressista.
  • congress: número da legislatura.
  • dwnom1: DW-NOMINATE score (dimensão 1).
  • dwnom2: DW-NOMINATE score (dimensão 2).
• Os dados vão da 80o legislatura (1947-48) ate a 112o legislatura (2011-2012)

• Vamos começar fazendo um gráfico de dispersão, scatter plot, para a 80o e 112o legislaturas. Republicanos e democratas serão caracterizados pela cor e pelo formato dos pontos.

plot_80 <- congress %>%
  filter(congress == 80) %>%
  ggplot(aes(x=dwnom1, y=dwnom2)) +
  geom_point(aes(shape=party, color = party), show.legend = FALSE) +
  scale_color_manual(values = c(Democrat = "blue",
                                Republican = "red",
                                Other = "green")) +
  scale_shape_manual(values = c(Democrat = "square",
                                Republican = "triangle",
                                Other = "circle")) +
  annotate("text",x=-1, y=1, label="Democrats", color = "blue") +
  annotate("text",x=0.5, y=-1, label="Republicans", color = "red") +
  scale_y_continuous("Racial liberalism/conservatism", limits = c(-1.5, 1.5)) +
  scale_x_continuous("Economic liberalism/conservatism", limits = c(-1.5, 1.5)) +
  ggtitle("80th Congress") + 
  coord_fixed() +
  theme_classic()

plot_112 <- congress %>%
  filter(congress == 112) %>%
  ggplot(aes(x=dwnom1, y=dwnom2)) +
  geom_point(aes(shape=party, color = party), show.legend = FALSE) +
  scale_color_manual(values = c(Democrat = "blue",
                                Republican = "red",
                                Other = "green")) +
  scale_shape_manual(values = c(Democrat = "square",
                                Republican = "triangle",
                                Other = "circle")) +
  annotate("text",x=-1, y=1, label="Democrats", color = "blue") +
  annotate("text",x=0.5, y=-1, label="Republicans", color = "red") +
  scale_y_continuous("Racial liberalism/conservatism", limits = c(-1.5, 1.5)) +
  scale_x_continuous("Economic liberalism/conservatism", limits = c(-1.5, 1.5)) +
  ggtitle("80th Congress") + 
  coord_fixed() +
  theme_classic()

grid.arrange(plot_80, plot_112, ncol=2)

• Note que as diferenças entre os partidos nas questões raciais são bem mais importantes no 80o Congresso do que no 112o Congresso. Vale o contrário nas questões econômicas.

• O diagrama de dispersão compara duas variáveis medidas para o mesmo conjunto de unidades em um gráfico com o valor de uma variável em um eixo e o de outra variável no outro eixo.

• O próximo exercício será comparar a mediana da primeira dimensão do DW-NOMINATE (liberalismo/conservadorismo econômico) em cada partido nas diversas legislaturas.
• Para isso vamos calcular as medianas e depois fazer um gráfico de linha.

• Calcular a mediana

median_dw1 <- congress %>%
  filter(party %in% c("Republican", "Democrat")) %>%
  group_by(party, congress) %>%
  summarize(median_dw1 = median(dwnom1)) %>%
  ungroup()

## `summarise()` has grouped output by 'party'. You can override using the
## `.groups` argument.

median_dw1

## # A tibble: 66 × 3
##    party    congress median_dw1
##    <chr>       <int>      <dbl>
##  1 Democrat       80     -0.126
##  2 Democrat       81     -0.207
##  3 Democrat       82     -0.179
##  4 Democrat       83     -0.174
##  5 Democrat       84     -0.223
##  6 Democrat       85     -0.231
##  7 Democrat       86     -0.259
##  8 Democrat       87     -0.258
##  9 Democrat       88     -0.285
## 10 Democrat       89     -0.292
## # ℹ 56 more rows

ggplot(data=median_dw1, aes(x=congress, y=median_dw1, color=party)) +
  geom_line() +
  labs(title = "DW-NOMINATE score (1st dimension)", x="Congress", y="Party") +
  theme_classic()

• O gráfico ilustra o distanciamento entre os dois partidos em termos de política econômica, os democratas ficam mais liberais (esquerda) e os republicanos mais conservadores (direita). A mudança dos republicanos é mais acentuada que a dos democratas.
• Esse fenômeno é chamado de polarização política.

• O que causou a polarização política? Essa é uma pergunta difícil de responder, o tema está agenda de pesquisa e gera debates na academia e outros meios. Alguns autores defendem que o aumento da desigualdade de renda está relacionado a essa polarização.
• Para medir a desigualdade de renda é comum usar o coeficiente de Gini.

• Para entender o coeficiente de Gini considere um gráfico onde o eixo horizontal representa o acumulado de pessoas em uma área específica (país, estado, etc) ordenado da menor para maior renda e o eixo vertical representa o acumulado de renda pertencente a indivíduos com renda igual ou menor que a da pessoa em um dado percentil de renda.
• A curva de Lorenz liga os pontos desse plano.
• Em uma distribuição de renda perfeitamente igualitária a curva de Lorenz é uma reta de 45o, isso ocorre porque x% da população detém, x% da renda da região. Essa reta é chamada de linha de igualdade.

• O coeficiente de Gini pode ser definido como a área entre a linha de igualdade e a curva de Lorenz dividida pela área abaixo da linha de igualdade. \[\small \mbox{coeficiente de Gini}=\frac{\mbox{área entre a linha de igualdade e a curva de Lorenz}}{\mbox{área abaixo da linha de igualdade}} \normalsize\]

• Na figura:

\[\mbox{coeficiente de Gini}=\frac{\mbox{área A}}{\mbox{área A + área B}}\]

• O coeficiente de Gini mede o grau de desigualdade de renda em uma sociedade. O coeficiente de Gini está entre zero (distribuição perfeitamente igualitária) e um (uma pessoa detém toda a renda).

• Para avaliar a relação entre desigualdade e polarização política comecemos comparando a evolução da medida de cada uma no tempo.
• Como medida de polarização usaremos a diferença entre as medianas dos dois partidos, como medida de desigualdade usaremos o coeficente de Gini.

polarization <- median_dw1 %>%
  pivot_wider(names_from = party,
              values_from = median_dw1) %>%
  mutate(polarization = Republican - Democrat)

polarization %>%
  ggplot(aes(congress, polarization)) +
  geom_point() +
  labs(title = "Politican polarization",
       x="Congress",
       y="Pepublican median -\n Democratic median") +
  theme_classic()

USGini <- read.csv("USGini.csv")
head(USGini)

##   X year  gini
## 1 1 1947 0.376
## 2 2 1948 0.371
## 3 3 1949 0.378
## 4 4 1950 0.379
## 5 5 1951 0.363
## 6 6 1952 0.368

USGini %>%
  ggplot(aes(x=year, y=gini)) +
  geom_point() +
  labs(title = "Gini coefficent",
       x= "Year",
       y = "Income inequality") +
  theme_classic

• É fácil ver que as duas variáveis crescem no tempo, mas lembre isso pode acontecer por diversos motivos. Nunca esqueça: correlação não é causalidade!
• Por falar em correlação…

• Correlação, ou coeficiente de correlação, é uma das medidas mais usadas para representar relações entre duas variáveis. A correlação representa como, em média, duas variáveis se movimentam conjuntamente em relação às próprias médias.
• Grosso modo, a correlação diz se quando uma variável está acima (abiaxo) da média a outra também está.

• Antes de uma definição formal de correlação é útil definir o valor-z (z-score). O valor-z representa o número de desvios-padrão uma observação está acima ou abaixo da média. \[\mbox{valor-z de x_i}=\frac{x_i-\mbox{média de x}}{\mbox{desvio padrão de x}}\]

• O valor-z padroniza a variável de forma que tranformações lineares não mudam o valor-z de uma variável, por exemplo, mudanças na unidade de medida não mudam o valor-z. \[\mbox{valor-z de }(ax_i+b)=\frac{(ax_i+b)-\mbox{média de }(ax_i+b)}{\mbox{desvio-padrão de }(ax_i+b)}\] \[=\frac{a\times(x_i-{\mbox{média de }x})}{a\times\mbox{desvio-padrão de }x}\] \[=\mbox{valor-z de }x\]

• O valor-z da i-ésima observação de uma variável \(x\) mede o número de desvios padrões que a observação está abaixo ou acima da média. É definido por: \[\mbox{valor-z de }x_i=\frac{x_i-\bar{x}}{S_x}\] Onde \(\bar{x}\) e \(S_x\) são, respectivamente, a média e o desvio-padrão de \(x\). O valor-z é uma medida de desvio da média e não é modificado por transformações lineares em \(x\).

• Conhecido o valor-z, a correlação é definida por: \[\mbox{correlação}(x,y)=\frac{1}{n}\sum_{i=1}^n(\mbox{valor-z de }x_i \times \mbox{valor-z de }y_i)\]
• Assim como no desvio padrão, é comum usar \(n-1\) no lugar de \(n\) no denominador, continua valendo que essa mudança tem pouco efieto quando \(n\) é grande.

• Correlação, ou coeficiente de correlação, mede o grau em duas variáveis estão associadas uma a outra. É definido como: \[\mbox{correlação}(x,y)=\frac{1}{n}\sum_{i=1}^n\left(\frac{x_i-\bar{x}}{S_x}\times\frac{y_i-\bar{y}}{S_y}\right)\mbox{ ou }\\\frac{1}{n-1}\sum_{i=1}^n\left(\frac{x_i-\bar{x}}{S_x}\times\frac{y_i-\bar{y}}{S_y}\right)\] Onde \(\bar{x}\) e \(\bar{y}\) são as médias e \(S_x\) e \(S_y\) são os desvios-padrão das variáveis \(x\) e \(y\), respectivamente.Correlação vai de -1 a 1 e não depende da unidade de medida da variável.

• No R a função cor() calcula a correlação.
• Para calcular a correlação entre polarização política e coeficiente de Gini precisamos primeiro arrumar os dados, pois o coeficiente de Gini engloba todos os anos e cada legislatura dura por dois anos. Serão usados os coeficientes de Gini dos anos pares, correspondem ao segundo ano de cada legislatura.

gini_2yr <- USGini %>%
  filter(row_number() %% 2 ==0) %>%
  select(gini) %>%
  pull()

pol_annual <- polarization %>%
  select(polarization) %>%
  pull()

cor(gini_2yr, pol_annual)

## [1] 0.9418128

• A correlação é alta e positiva. Isso sugere que desigualdade de renda e polarização política andam juntas.
• Podemos falar de causalidade? Não, muitas variáveis mostram tendência de crescimento nesse período e, por óbvio, nem todas possuem relação de causalidade.

• Algumas vezes estamos interessados em comparar toda a distribuição de duas variáveis, não apenas medidas como média ou mediana, uma alternativa para fazer essa comparação é por meio do histograma de cada distribuição.
• Nesses tipos de comparação é importante que as escalas dos gráficos sejam as mesmas
• A função facet_grid() do ggplot2 permite fazer gráficos do tipo. No exemplo vamos comparar as distribuições da segunda dimensão do DW-NOMINATE score entre republicanos e democratas da 112o legislatura.

congress %>%
  filter(congress == 112, party %in% c("Republican", "Democrat")) %>%
  ggplot(aes(x=dwnom2, y=after_stat(density))) +
           geom_histogram(binwidth = .2) +
           facet_grid(party ~ .) +
           labs(x="Racial liberalism/conservatism dimension",
                y="Density") +
           theme_classic()

• As distribuições são semelhantes, mas é possível notar que a distribuição dos democratas tem uma cauda superior mais longa do que a dos republicanos (mais observações com valores mais altos).
• Outra diferença é que a distribuição dos republicanos está mais concentrada no centro do que a dos democratas.

• O gráfico quantil-quantil, ou gráfico Q-Q, é outra forma de comparar distribuições.
• O gráfico Q-Q é um diagrama de dispersão onde os pontos representam os mesmo quantis das distribuições, por exemplo, um ponto referente a mediana de uma distribuição também é referente a mediana da outra distribuição.

• Se duas distribuiições são idênticas, então os quantis serãoos mesmos e o gráfico Q-Q será uma reta de 45 graus.Pontos acima da reta de 45 graus signficam que o valor daquele quantil na variável do eixo vertical é maior do que o do variável no eixo horizontal. Vale o contrário para pontos abaixo da reta de 45 graus.
• A inclinação da curva nos dá ideia da dispersão das duas variáveis. Se a curva for mais inclinada do que a reta de 45 graus, então a distribuição do exio vertical tem mais dispersão do que a do eixo horizontal, vale o contrário quando a inlcinação da curva é menor do que 45 graus.

• O script seguinte mostra a construção do gráfico QQ para as distribuições da segunda dimensão do DW-NOMINATE score entre republicanos e democratas.
• Repare que o exemplo usa a função geom_point() para criar o gráfico.

quantile_probs <- seq(from=0, to=1, by=0.01)
quantile_names <- as.character(quantile_probs)

quantiles <- congress %>%
  filter(congress == 112) %>%
  group_by(party) %>%
  summarize(dwnom_quantile = quantile(dwnom2, probs=quantile_probs),
            quantile = quantile_names) %>%
  pivot_wider(names_from=party, values_from=dwnom_quantile) %>%
  ungroup()

quantiles %>%
  ggplot(aes(x=Democrat, y=Republican)) +
  geom_point(shape=1) +
  ylim(-1.5, 1.5) +
  xlim(-1.5, 1.5) +
  geom_abline(intercept=0, slope=1) +
  ggtitle("Racial liberalism/conservadorism dimension") +
  coord_fixed()

• A função qqplot() do R básico oferece uma maneira mais simples de fazer um gráfico Q-Q, porém sem as comodidades do tidyverse

dem112 <- filter(congress, party == "Democrat", congress == 112)
rep112 <- filter(congress, party == "Republican", congress == 112)

qqplot(x=dem112$dwnom2,
       y=rep112$dwnom2,
       xlab = "Democrats",
       ylab = "Republicans",
       xlim = c(-1.5, 1.5), ylim = c(-1.5, 1.5),
       main = "Racial liberalism/conservadorism dimension")

• Comno os pontos representando os quantis menores estão acima da linha de 45 graus podemos dizer que os liberais republicanos são mais conservadores do que os liberais democratas. Isso é verdade porque esses quantis tem maiores valores (ou seja, mais conservadores) para os republicanos do que os quantis correspondetes para os democratas.
• Por outro lado, os pontos correspondentes aos maiores quantis estão abaixo da linha de 45 graus, ou seja, nos quantis maiores, os mais conservadores, os democratas tem valores maiores (mais conservadores) do que os republicanos. Os conservadores democratas são mais conservadores do que os conservadores republicanos.

• Vimos que na 1120 legislatura, Democratas e Republicanos formaram dois grupos ideologicamente distintos, mas será que existem agrupamentos (clusters) ideológicos dentro do mesmo partido?
• Como encontrar esses grupos?
• Antes de responder essas perguntas é útil falar sobre matrizes no R

• Assim como os objetos da classe data.frame, os objetos da classe matrix são estruturas retângulares. Os dois tipos de objeto possuem várias semelhanças, mas apresentam algumas diferenças relevantes.
• A mais importante, enquanto data.frame admite diferentes tipos de variáveis (numéricas, textos, lógicas, etc), o objeto da classe matrix só admite um tipo de variável, normalmente numérica.
• Matrizes são particularmente úteis para aplicações em álgebra linear e são computaiconalmente mais eficientes do que data.frame.

• Um objeto da classe matrix é criado pela função matrix(), para isso é necessário fornecer um vetor com os valores dos elementos e informar se o vetor deve ser lido por linhas ou por colunas bem como o número de linhas/colunas.
• É possível dar nomes as linhas e/ou colunas com as funções rownames() e colnames().

x <- matrix(1:12, nrow=3, ncol=4, byrow=TRUE)
dim(x)

## [1] 3 4

x

##      [,1] [,2] [,3] [,4]
## [1,]    1    2    3    4
## [2,]    5    6    7    8
## [3,]    9   10   11   12

x <- matrix(1:12, nrow=3, ncol=4, byrow=TRUE)
rownames(x) <- c("a", "b", "C")
colnames(x) <- c("d", "e", "f", "g")
x

##   d  e  f  g
## a 1  2  3  4
## b 5  6  7  8
## C 9 10 11 12

• Como, ao contrário do data.frame as matrizes só admitem uma classe de variável, se alguém trasformar uma data.frame em matriz com a função as.matrix() o R vai transformar a classe das variáveis de forma que todas fiquem com a mesma classe.

y <- data.frame(y1 = as.factor(c("a", "b", "c")),
                y2 = c(0.1,0.2,0.3))
str(y)

## 'data.frame':    3 obs. of  2 variables:
##  $ y1: Factor w/ 3 levels "a","b","c": 1 2 3
##  $ y2: num  0.1 0.2 0.3

z <- as.matrix(y)
z

##      y1  y2   
## [1,] "a" "0.1"
## [2,] "b" "0.2"
## [3,] "c" "0.3"

• Algumas funções permitem realizar operações com as linhas ou colunas de uma matriz, por exemplo, colSums() e rowMeans().

colSums(x)

##  d  e  f  g 
## 15 18 21 24

rowMeans(x)

##    a    b    C 
##  2.5  6.5 10.5

• Outra importante classe de objetos no R são as listas. Listas podem aramazenar diferentes tipos de objetos e, ao contrário do data.frame, com diferentes dimensões, por exemplo, uma lista pode armazenar um vetor de caracteres e outro de números com tamanhos diferentes. Listas podem armazenar até mesmo data.frames de tamanhos diferentes.
• Elementos de uma lista podem ser acessados com o síbolo $, como variáveis em um data.frame ou com [[]] tendo como argumento o número ou nome do objeto que desejamos acessar.

x <- list(y1 = 1:10, y2 = c("hi", "hello", "hey"),
          y3 = data.frame(z1=1:3, z2 = c("good", "bad", "ugly")))

x$y1

##  [1]  1  2  3  4  5  6  7  8  9 10

x[[2]]

## [1] "hi"    "hello" "hey"

x[["y3"]]

##   z1   z2
## 1  1 good
## 2  2  bad
## 3  3 ugly

• A função names() se aplicada em uma lista retorna os nomes dos objetos na lista e a função length() retorna o número de objetos na lista.

names(x)

## [1] "y1" "y2" "y3"

length(x)

## [1] 3

• O algoritmo k-means é uma forma popular de agrupar dados. O objetico do algoritmo é distribuir as variáveis em k grupos semelhantes onde cada grupo é associado a um centro. O primeiro passo é definir o número de centros e “espalhar” os centros no espaço de dados, depois os dados são classficados pela distância desse centro. Após a classificação dos dados, novos centros são calculados e os dados são reclassficados de acordo com os novos centros. O algoritmo para quando nenhum dados mudar de classificação de uma rodada para outra.

• No link https://www.tidymodels.org/learn/statistics/k-means/ você encontra uma animação que ilustra bem o algoritmo k-means.

• O algoritmo k-means produz um determinado número de grupos, ou clusters, chame esse número de \(k\), e consiste nos seguintes passos:
  • Escolha os \(k\) centros iniciais;
  • Dados os centros, classifique cada observação no cluster com centro mais próximo da observação;
  • Escolha os novos centros de cada cluster de forma que o centro seja a média de cada coordenada das observações pertencentes ao cluster;
  • Repita os dois passos anteriores até que nenhuma observação mude de grupo quando os centros são calculados.

• É costume padronizar os dados antes de aplicar o algoritimo k-means, isso é feito usando o z-valor, ou seja, dimunidno da média e dividindo pelo desvio padrão.
• A padronização dos dados elimina problemas de escala.

• Como exemplo vamos aplicar o algortimo k-means para as 80o e 120o legislaturas americanas usando \(k=2\) e \(k=4\).
• A função kmeans() aplica o algoritmo, entre seus argumentos estão centers que determina o número de clusters, iter.max que é o número máximo de iterações e nstart o número de centros aleatórios que devem ser escolhidos.
• É recomendado especificar o argumento nstart para que o algoritimo seja aplicado diversas vezes com diferentes pontos de partida, a função kmeans() retorna o melhor agrupamento.

k80two.out <- congress %>%
  filter(congress == 80) %>%
  select(dwnom1, dwnom2) %>%
  kmeans(centers = 2, nstart = 5)

k112two.out <- congress %>%
  filter(congress == 112) %>%
  select(dwnom1, dwnom2) %>%
  kmeans(centers = 2, nstart = 5)

• Os objetos k80two.out e k112two.out são listas contendo várias informações retornadas pela função kmeans(), entre as quais o número de iterações, iter, um vetor dizendo a que cluster pertence cada observação, cluster e uma matriz com os centros, centers.

names(k80two.out)

## [1] "cluster"      "centers"      "totss"        "withinss"     "tot.withinss"
## [6] "betweenss"    "size"         "iter"         "ifault"

• Nos dois casos o algortimo convergiu com apenas uma iteração, é importante checar esse número pois o default para o número máximo de iterações é dez, se o algortimo usou as dez iterações é necessário tentar de novo com o argumento iter definido como um número maior do que dez.

k80two.out$iter

## [1] 1

k112two.out$iter

## [1] 1

• O próximo passo é mostrar os centros de cada cluster.

k80two.out$center

##        dwnom1     dwnom2
## 1  0.14681029 -0.3389293
## 2 -0.04843704  0.7827259

k112two.out$center

##       dwnom1     dwnom2
## 1  0.6776736 0.09061157
## 2 -0.3912687 0.03260696

• Para seguir com o exemplo é útil usar a função tidy() do pacote tidymodels ou broom, para transformar os resultados em data.frames fáceis de manipular.

library(tidymodels)

k80two.clusters <- tidy(k80two.out)
k80two.clusters

## # A tibble: 2 × 5
##    dwnom1 dwnom2  size withinss cluster
##     <dbl>  <dbl> <int>    <dbl> <fct>  
## 1  0.147  -0.339   311     54.9 1      
## 2 -0.0484  0.783   135     10.9 2

k112two.clusters <- tidy(k112two.out)
k112two.clusters

## # A tibble: 2 × 5
##   dwnom1 dwnom2  size withinss cluster
##    <dbl>  <dbl> <int>    <dbl> <fct>  
## 1  0.678 0.0906   242     27.1 1      
## 2 -0.391 0.0326   201     38.8 2

• Quantos Democratas e quantos Republicanos em cada cluster na 80o legislatura?

congress %>%
  filter(congress == 80) %>%
  mutate(cluster2 = k80two.out$cluster) %>%
  group_by(party, cluster2) %>%
  count() %>%
  pivot_wider(names_from = cluster2, values_from = n)

## # A tibble: 3 × 3
## # Groups:   party [3]
##   party        `1`   `2`
##   <chr>      <int> <int>
## 1 Democrat      62   132
## 2 Other          2    NA
## 3 Republican   247     3

• Quantos Democratas e quantos Republicanos em cada cluster na 112o legislatura?

congress %>%
  filter(congress == 112) %>%
  mutate(cluster2 = k112two.out$cluster) %>%
  group_by(party, cluster2) %>%
  count() %>%
  pivot_wider(names_from = cluster2, values_from = n)

## # A tibble: 2 × 3
## # Groups:   party [2]
##   party        `2`   `1`
##   <chr>      <int> <int>
## 1 Democrat     200    NA
## 2 Republican     1   242

• Para o próximo exercício serão usados quatro clusters e o resultado será ilustrado em gráficos.

k80four.out <- congress %>%
  filter(congress == 80) %>%
  select(dwnom1, dwnom2) %>%
  kmeans(centers = 4, nstart = 5)

k112four.out <- congress %>%
  filter(congress == 112) %>%
  select(dwnom1, dwnom2) %>%
  kmeans(centers = 4, nstart = 5)

k80four.cluster <- tidy(k80four.out)
k112four.cluster <- tidy(k112four.out)

congress80 <- filter(congress, congress == 80) %>%
  mutate(cluster4 = factor(k80four.out$cluster))
congress112 <- filter(congress, congress == 112) %>%
  mutate(cluster4 = factor(k112four.out$cluster))

ggplot() +
  geom_point(data = congress80, aes(dwnom1, dwnom2, color=cluster4)) +
  geom_point(data = k80four.cluster, aes(dwnom1, dwnom2), size=3, shape=8) +
  ylim(-1.5,1.5) +
  xlim(-1.5,1.5) +
  coord_fixed() +
  theme(legend.position = "none")

ggplot() +
  geom_point(data = congress112, aes(dwnom1, dwnom2, color=cluster4)) +
  geom_point(data = k112four.cluster, aes(dwnom1, dwnom2), size=3, shape=8) +
  ylim(-1.5,1.5) +
  xlim(-1.5,1.5) +
  coord_fixed() +
  theme(legend.position = "none")

• A análise com cluster sugere que a questão racial define clusters no partido e a questão econômica separa os partidos.
• Na 80o legislatura fica bem clara a divisão entre os democratas.

• Algoritmos de agrupamento, como o k-means, são exemplos de métodos de aprendizado não supervisionado. Ao contrário dos métodos supervisionados não há uma variável de resultado especificada com antecedência. De fato, o objetivo do aprendizado não supervisionado é descobrir estrutras não antecipadas dos dados, como, por exemplo, clusters.

• A dificuldade do aprendizado não supervisionado é não ter uma medida clara de sucesso ou fracasso, é difcíl saber se o algoritmo está produzindo o resultado “correto”. Por exemplo, não sabemos se quatro é o número certo de clusters no nosso exemplo.
• Por esta razão, o julgamento humano é fundamental para avaliar o resultado de um algoritmo de agrupamento.