© Ricardo Solar/UFMG - compartilhamento e utilização não-comercial livres. Não reproduzir sem autorização >http://doi.org/10.5281/zenodo.7392285

Agora, começamremos a fazer as análises de maneira mais direta no R. Vamos utilizar uma tabela que contém dados de altura e peso de indivíduos. Os dados foram coletados em turmas de R em vários anos, de turmas de Graduação e Pós-Graduação. Vamos primeiro importar os dados, detalhando como isso deve ser feito.

1 Importação de dados de uma planilha preexistente

Ainda que saibamos que tudo que faremos neste curso é importante, essa parte merece uma atenção especial! Se você não conseguir importar seus dados de maneira correta e direta ao R, de muito pouca utilidade o programa será pra você! O carregamento de dados pode parecer demasiado técnico do ponto de vista computacional, e por isso farei algumas recomendações que podem facilitar esse processo. Aqui é importante destacar que tais recomendações não são as únicas formas de carregar seus dados, mas padronizar de alguma maneira essa tarefa pode diminuir suas dores de cabeça.

De maneira simplificada, importar seus dados no R envolve a seguinte sequência:

  1. criar ou baixar sua planilha de dados. Essa planilha terá, como padrão geral, as observações nas linhas e a as variáveis nas colunas, de acordo com o exemplo abaixo:
ID Var_1 Var_2 Var 3
amostra_1 valor_v1 valor_v2 valor_v3
amostra_2 valor_v1 valor_v2 valor_v3
amostra_n valor_v1 valor_v2 valor_v3
  1. localizar e avisar o R em qual pasta o arquivo de trabalho está localizado - funções getwd e setwd que veremos abaixo.
  2. importar os dados no R - funções de importação de dados (read.table, read.csv, read.excel, etc)
  3. conferir se os dados foram importados corretamente, de acordo com o seu conhecimento dos dados - observação do Environment e funções summary, str, plot, View, etc.

1.1 Formatos de importação de dados

Os dados podem ser importados para o R em vários formatos, basicamente todas as formas de texto existentes. Exatamente por esse motivo, recomendo que você escolha um formato preferencial que use sempre, evitando se confundir. É importante, todavia, ter uma noção dos principais formatos de dados que trabalha-se no R. Vamos ver alguns dos mais comuns:

  • .txt - o formato TXT é bastante versátil. É um formato de texto simples, onde os valores podem ser separados por espaços ou tabulações. É um arquivo de texto puro, isto é, sem formatação ou fórmulas. Isso é muito útil, pois esse tipo de arquivo não será importado para o R caso haja qualquer erro remanescente. Para uso no R, recomendo fortemente exportar o arquivo do Excel (ou seus equivalentes) como texto (separado por tabulações). A importação desse tipo de arquivo no R é bem direta, através da função read.table.
  • .csv - o formato é uma abreviação para comma separated values (de valores separados por vírgula). Também é um arquivo de texto puro, como o TXT. Há alguns detalhes, se arquivo CSV for exportado de uma versão que utiliza formatação em inglês do Excel, os valores em cada linha serão separados por vírgulas (,). Se o arquivo CSV for gerado usando uma versão que usa formatação dinamarquesa do Excel, a vírgula é usada para casas decimais e os valores são separados por ponto e vírgula (;). Essa é uma das questões técnicas que mencionei e que podem enrolar o uso de alguns arquivos. A importação desse tipo de arquivo no R é bem direta, através da função read.csv ou read.csv2.
  • xls ou xlsx - é possível ainda importar os dados diretamente de um arquivo Excel. O R não faz essa importação direta e você precisará lançar mão de um pacote (veremos em breve) para importar esses arquivos (por exemplo o pacote readxl, com sua função read.excel). A importação direta do arquivo Excel é bem interessante, mas o arquivo vem carregado de formatações, fórmulas e adjuntos trazidos pelo Excel. É importante tomar cuidado na conferência dos dados.

1.1.1 O que observar antes de importar seus dados no R?

Há informações que você precisa saber antes de importar qualquer dado no R, que baseiam-se majoritariamente em como a tabela foi construída, vou listar alguns abaixo, mas vale a pena ir construindo seu formato de fazer a tabela para evitar problemas.

  • A primeira linha contém nomes de variáveis (o chamado cabeçalho)?
  • Há valores ausentes? Eles estão listados de maneira padronizada como NA, tanto no formato, quanto maiúsculas?
  • Os números com casas decimais são listados usando pontos (.) ou vírgulas (,)? Todas as casas decimais utilizam o mesmo separador? Aqui mora um dos maiores problemas de digitação de dados em R. O padrão do programa é separar decimais por pontos, então se usou pontos, não precisa mais se preocupar.
  • caso haja dados de variáveis que devam ser lidos como Factors? Se houver, é necessário declarar isso na importação como stringsAsFactors = T, uma vez que o padrão do programa é que eles sejam lidos como texto.
  • O formato de arquivo é TXT, CSV ou Excel? No caso de estar usando CSV, é separado por vírgula ou ponto-e-vírgula? Sei trabalhar com o formato ou vale a pena abrir no Excel e exportar em outro formato?

Agora que sabemos tudo isso, e estamos certos do nosso arquivo, posso continuar a importar meus dados.

1.2 Localizar e mudar o diretório de trabalho

Para trabalhar com os dados no R, precisamos informar qual é o nosso diretório de trabalho. A forma mais simples de fazer isso, como já vimos, é criar uma pasta para cada análise que faremos e trabalhar dentro daquela pasta. Essa prática apresenta uma série de vantagens. A maior delas é saber que todos os arquivos que preciso importar estão na mesma pasta. Outra grande vantagem é que todo arquivo que eu exportar, será gravado naquela pasta também.

1.2.1 Mudando o diretório - setwd

Nessa etapa, é importante ter alguma intimidade com seu computador. Isso vai fazer sua vida mais fácil.

Há alguns caminhos possíveis, mas todos eles devem em última instância lhe gerar um caminho no computador que possa ser colocado no script, para que em rodadas futuras daquela análise isso já esteja resolvido.

  1. Clicar no Menu e Session > Set Working Directory > Choose Directory. Isso abrirá uma janela para que você navegue nos seus arquivos até a pasta onde estão suas planilhas. Veja que você escolherá a pasta, não o arquivo.
  2. No painel Files, navegue até a pasta onde estão seus arquivos. Ao chegar na pasta, sem clicar no arquivo, há um ícone no canto superior direito do painel com uma engrenagem azul escrito More. Clique em More > Set As Working Directory.

Estes dois comandos irão colocar num console um comando começado com setwd similar ao que apareceu aqui no endereço do meu computador: setwd("~/Dropbox/UFMG/Disciplinas/R UFAC/Aulas/teste t e anova"). Veja que este endereço é exclusivo desta pasta no meu computador, haverá um diferente no seu. Copie este endereço do console e cole ele no seu script, para que em sua próxima rodada dessa análise não haja necessidade de fazer isso de novo.

  1. Uma última forma é, caso já possua conhecimento da organização de dados no seu computador, digite esse endereço diretamente no script, isso vai salvar muito tempo evitando cliques desnecessários.

1.3 Mãos à obra!

Vamos agora importar nosso arquivo, que será utilizado nas análises a seguir. Conferindo:

  • Está em formato de texto, separado por tabulações. *.txt;
  • A primeira linha contém nomes de variáveis (o chamado cabeçalho);
  • Não há valores ausentes, portanto não preciso me preocupar com NAs
  • Os números com casas decimais são listados usando pontos, e todas as casas decimais utilizam o mesmo separador. Vou conferir isso novamente depois.
  • Tenho dados que precisam ser entendidos como Factors

Com isso em mente, vamos lá.

Primeira parte, mudar o diretório:

setwd("~/Dropbox/UFMG/Disciplinas/R UFAC/Aulas/teste t e anova")

##caso eu queira conferir, posso usar os comandos
getwd() #para que o programa me mostre o diretório que ele está setado
## [1] "/Users/ricardosolar/Dropbox/UFMG/Disciplinas/R UFAC/Aulas/teste t e anova"
##e também
dir() ##para me mostrar os arquivos que estão na pasta que estou trabalhando.
##  [1] "dados_qualif.txt"                      
##  [2] "Herbivoria.txt"                        
##  [3] "Regressão Linear simples e ANCOVA.Rmd" 
##  [4] "Regressão-Linear-simples-e-ANCOVA.html"
##  [5] "Regressão-Linear-simples-e-ANCOVA.log" 
##  [6] "Regressão-Linear-simples-e-ANCOVA.tex" 
##  [7] "rsconnect"                             
##  [8] "RTurmas.txt"                           
##  [9] "Teste T e F - rotina usual.Rmd"        
## [10] "Teste-T-e-F---rotina-urual.html"       
## [11] "Teste-T-e-F---rotina-usual.html"

Agora posso importar os dados, que como vimos, estão com o nome RTurmas.txt. Como os dados estão em txt, vou usar o comando read.table, e detalhar os passos que usamos na função abaixo:

read.table( #comando base
  "RTurmas.txt", #arquivo a ser importado
  header = TRUE, #Informa que há um cabeçalho, ou seja, a primeira linha deve ser considerada como nomes dos dados, e não deve ser entendida como parte em si dos dados
  stringsAsFactors = TRUE, #informo que, caso haja texto que se repita nas variáveis, devem ser considerados como fatores
  dec = "." #como havia dito, no nosso caso esse comando é dispensável, mas deixo aqui registrado que, caso sua decimal seja ",", é aqui que deve ser informada.
  )

O problema desse comando anterior é que ele só leu a planilha!! Para armazenar e utilizar no R, vale lembrar das lições anteriores … devemos designá-lo a um objeto, é bem simples, e fica assim no uso corriqueiro

dados <- read.table("RTurmas.txt", header = TRUE, stringsAsFactors = TRUE)

E pra poder conferir os dados, temos várias opções, vamos explorar algumas:

  • summary(x) - dá um resumo geral dos dados, focando nos valores e métricas das variáveis
summary(dados)
##      ID_Ind          Altura           Peso        Sexo    Turma   
##  Min.   : 1.00   Min.   :150.0   Min.   : 44.50   F:45   Grad:45  
##  1st Qu.:23.75   1st Qu.:163.0   1st Qu.: 56.77   M:47   Pos :47  
##  Median :46.50   Median :168.0   Median : 66.02                   
##  Mean   :46.50   Mean   :170.1   Mean   : 68.33                   
##  3rd Qu.:69.25   3rd Qu.:177.2   3rd Qu.: 76.03                   
##  Max.   :92.00   Max.   :195.0   Max.   :113.70                   
##       Ano      
##  Min.   :2015  
##  1st Qu.:2016  
##  Median :2017  
##  Mean   :2017  
##  3rd Qu.:2017  
##  Max.   :2019
  • str(x) - resume a estrutura das variáveis
str(dados)
## 'data.frame':    92 obs. of  6 variables:
##  $ ID_Ind: int  1 2 3 4 5 6 7 8 9 10 ...
##  $ Altura: int  162 160 165 165 164 165 158 160 166 165 ...
##  $ Peso  : num  44.5 46 52 53 54.2 56.1 57.6 58 62 79.1 ...
##  $ Sexo  : Factor w/ 2 levels "F","M": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Turma : Factor w/ 2 levels "Grad","Pos": 1 1 1 1 2 2 2 2 2 2 ...
##  $ Ano   : int  2015 2015 2015 2015 2015 2015 2015 2015 2015 2015 ...
  • dim(x) - mostra os números de linhas e colunas dos dados
dim(dados)
## [1] 92  6
  • names(x) ou colnames(x) - mostra os nomes das colunas dos dados, muito útil para ver se isso foi importado corretamente
names(dados)
## [1] "ID_Ind" "Altura" "Peso"   "Sexo"   "Turma"  "Ano"
colnames(dados)
## [1] "ID_Ind" "Altura" "Peso"   "Sexo"   "Turma"  "Ano"
  • head(x) e tail(x) - mostra o início (head) ou o final (tail) dos dados, como padrão tendo as 6 primeiras ou 6 últimas linhas.
#"cabeça" dos dados
head(dados)
#"cauda" dos dados
tail(dados)
  • plot(x) - por fim, posso fazer um plot geral dos dados, pra entender como estão distribuídos e, visualmente, identificar algum problema
plot(dados)

Com os dados conferidos e corretos, posso começar a trabalhar os mesmos.

2 Análise dos dados

2.1 Teste T simples

2.1.1 Pergunta a ser respondida

Queremos inicialmente verificar se a altura dos estudantes que cursaram a disciplina diferem entre os sexos biológicos. Isso pode ser feito através de um teste T.

A primeira premissa é que a variável resposta (Altura) deve ter distibuição normal. Como verificar isso? Uma das formas é o teste de Shapiro, que é implementado no R através de shapiro.test(x). Quero testar a variável Altura, então preciso informar ao programa naquele formato dados$Altura. O teste de shapiro tem como hipótese nula que a variável em questão possui distribuição Normal.

shapiro.test(dados$Altura)
## 
##  Shapiro-Wilk normality test
## 
## data:  dados$Altura
## W = 0.97605, p-value = 0.08814

O teste de shapiro não foi significativo, ou seja, não rejeito a hipótese nula. Portanto, minha variável possui distribuição Normal. Posso prosseguir com o teste T. Como vimos antes, o teste T é implementado através do comando t.test(y~x). Que irei armazenar num objeto chamado testeT.

(testeT <- t.test(Altura~Sexo, data=dados)) ## colocar todo o comando entre parêntesis permite que o R crie o objeto E mostre seu resultado
## 
##  Welch Two Sample t-test
## 
## data:  Altura by Sexo
## t = -9.6706, df = 89.763, p-value = 1.418e-15
## alternative hypothesis: true difference in means between group F and group M is not equal to 0
## 95 percent confidence interval:
##  -16.03837 -10.57155
## sample estimates:
## mean in group F mean in group M 
##        163.3333        176.6383

Veja que terminei o comando com data=dados. Fiz isso para informar ao programa que os termos Altura e Sexo estão dentro do dataframe dados.

Quando fiz o teste, gravei seus resultados no objeto testeT para poder extair dele os resultados sem ter que copiar e colar do console. Isso será principalmente útil para os gráficos. Vamos começar bem simples, fazendo um gráfico contendo o valor das médias e o intervalo de confiança. Para isso, já utilizarei o pacote GGPLOT2. Ainda não explicarei a ideia de semântica de gráficos, que virá logo em breve.

Vamos ver o que está contido dentro do objeto - usaremos o comando str(x)

str(testeT)
## List of 10
##  $ statistic  : Named num -9.67
##   ..- attr(*, "names")= chr "t"
##  $ parameter  : Named num 89.8
##   ..- attr(*, "names")= chr "df"
##  $ p.value    : num 1.42e-15
##  $ conf.int   : num [1:2] -16 -10.6
##   ..- attr(*, "conf.level")= num 0.95
##  $ estimate   : Named num [1:2] 163 177
##   ..- attr(*, "names")= chr [1:2] "mean in group F" "mean in group M"
##  $ null.value : Named num 0
##   ..- attr(*, "names")= chr "difference in means between group F and group M"
##  $ stderr     : num 1.38
##  $ alternative: chr "two.sided"
##  $ method     : chr "Welch Two Sample t-test"
##  $ data.name  : chr "Altura by Sexo"
##  - attr(*, "class")= chr "htest"

Ok, temos que o que nos interessa nessa estrtura são os pontos $estimate, onde estão armazenadas as médias e $stderr, que intuitivamente sabemos que armazena o erro padrão. Para fazer o gráfico, preciso montar isso num dataframe que permita ao programa entender os valores da forma correta. É bem simples, vamos lá:

dados.plot <- data.frame(medias = testeT$estimate, intervalo=abs(testeT$stderr), sexo=c("F", "M")) ##este último eu deduzi pela ordem alfabética, que o R sempre segue como padrão, F vem antes de M.

Veja acima que pedi a criação de um objeto data.frame, com o nome de dados.plot. Nele coloquei os elementos contidos em testeT$estimate e os elementos em testeT$stderr, além de criar uma outra coluna que continha os nomes c("F", "M"). Para incluir os pontos reais que deram origem às médias, devo usar o geom_jitter, mas preciso informar quque os dados estão em outro dataframe, confira abaixo. Veja que para o intervalo de confiança, me interessam somente os valores em módulo, ou valores absolutos.

library(ggplot2)
ggplot(dados.plot, aes(x=sexo, y=medias))+
  geom_jitter(data=dados, aes(x=Sexo, y=Altura), pch=18)+
  geom_point(pch=21, size=4, fill="grey50")+
  geom_errorbar(aes(ymin=medias-intervalo, ymax=medias+intervalo), position = position_dodge(), width=.3)+
  xlab("Sexo biológico")+
  ylab("Altura (cm)")

Assim então temos o teste T simples, sempre entre uma variável resposta e uma variável explicativa categórica (duas categorias).

2.2 Teste T Pareado

Para o teste T pareado, nosso interesse é diferente. Queremos saber **como a mesma amostra se comporta em dois momentos distintos, ou seja, quero acompanhar essas variáveis em dois momentos.

Para exemplificar o uso desse tipo de teste, vamos trabalhar com dados de um exemplo em que as informações foram coletadas antes e depois de plantas serem expostas à herbivoria. A hipótese é que a quantidade de compostos secundários das folhas aumenta após a exposição à herbivoria.

dadosherb <- read.table("Herbivoria.txt", h=T, stringsAsFactors = T)
summary(dadosherb)
##      Planta        Tempo     Concentracao  
##  Min.   : 1.0   after :10   Min.   : 8.70  
##  1st Qu.: 3.0   before:10   1st Qu.:11.97  
##  Median : 5.5               Median :13.45  
##  Mean   : 5.5               Mean   :13.06  
##  3rd Qu.: 8.0               3rd Qu.:14.30  
##  Max.   :10.0               Max.   :16.50

A nossa pergunta então é bem importante para definir o modelo. Meu interesse aqui é saber se a concentração de compostos aumenta depois do tratamento. Portanto, nosso modelo é a \(Concentração ~ Tempo\), levando em conta que eu já espero que os valores sejam maiores depois e que há valores pareados para a mesma planta antes e depois.

shapiro.test(dadosherb$Concentracao)
## 
##  Shapiro-Wilk normality test
## 
## data:  dadosherb$Concentracao
## W = 0.96403, p-value = 0.6271
(testeTpar <- t.test(Concentracao~Tempo, data=dadosherb, paired=TRUE, alternative="greater"))
## 
##  Paired t-test
## 
## data:  Concentracao by Tempo
## t = 2.272, df = 9, p-value = 0.0246
## alternative hypothesis: true mean difference is greater than 0
## 95 percent confidence interval:
##  0.1043169       Inf
## sample estimates:
## mean difference 
##            0.54

Nossos resultados, primeiro do teste de Shapiro, mostram que posso usar o teste T, pois a distribuição da variável é Normal. O restultado do teste T pareado, com a hipótese alternativa de que depois do tratamento os valores são maiores, é significativa (p=0.025).

Os gráficos de testes T pareados são um pouco diferentes, pois nos interessa ver o comportamento de cada ponto. Para isso, vamos usar um pacote que permite essa plotagem, tendo como base o pacote ggplot2.

Vou aproveitar esse momento para mostrar como filtar dados, usando o pacote library(dplyr). Observem que para fazer os gráfico, preciso separar os dados, no caso entre before e after. O comando filter permite filtrar um dataframe de acordo com uma condição determinada. No nosso caso, observe que filtrei before e depois after. Como só era do meu interesse a coluna `$Concentracao’, adicionei essa indexação ao final do comando.

# Plot paired data
library(PairedData)
before <- dplyr::filter(dadosherb, Tempo=="before")$Concentracao
after <- dplyr::filter(dadosherb, Tempo=="after")$Concentracao

pd <- paired(before, after)

plot(pd, type = "prof")+
  geom_line(col="coral", linewidth=1)+
  #geom_jitter(width=.1)
  geom_point(col="#801FA4")

3 Teste F e ANOVA

E se nosso interesse fosse testar mais de uma variável, ou variáveis com mais de dois níveis? Nessa situação, podemos lançar mão da Análise de Variância(ANOVA), com o teste F. De fato, veremos que é sempre possível usar esse teste, pois quando entrarmos em GLM, esse tipo de análise permite acomodar modelos não-normais, o que nos ajuda muito!

Uma coisa importante aqui é que inicialmente também temos a premissa de normalidade das populações de onde tiramos os dados. Poderíamos também seguir com o teste de Shapiro para checar a normalidade, mas temos algumas observações a fazer. Em geral, uma ANOVA de um teste bastante robusto contra violações da premissa de normalidade, desde que os tamanhos de amostra sejam suficientemente grandes.

Além disso, se você tiver tamanhos de amostra extremamente grandes, testes estatísticos como o teste de Shapiro-Wilk quase sempre dirão que seus dados não são normais. Por esse motivo, geralmente é melhor inspecionar seus dados visualmente, através de análise dos resíduos. Simplesmente olhando para os gráficos, você pode ter uma boa ideia se os dados são ou não normalmente distribuídos. Faremos isso consistentemente de agora em diante.

Agora que já sabemos isso tudo, podemos testar nossa nova pergunta de interesse com aquele primeiro conjunto de dados, que contém as informações das turmas de R.

Nossa pergunta agora é se o Peso da turma é função do Sexo biológico e do fato de estas pessoas estarem na Graduação ou na Pós. Nosso modelo agora passa a ser então Peso~Sexo+Turma. Para isso, usaremos o comando aov, e colocaremos esse modelo sob prova. Para testar as premissas do modelo, usaremos a análise de resíduos, implementada através do pacote library(DHARMa).

modelo <- aov(Peso~Sexo+Turma, data=dados)

library(DHARMa)
simulateResiduals(modelo, plot=T)

## Object of Class DHARMa with simulated residuals based on 250 simulations with refit = FALSE . See ?DHARMa::simulateResiduals for help. 
##  
## Scaled residual values: 0.204 0.252 0.48 0.524 0.132 0.152 0.2 0.212 0.356 0.904 0.34 0.448 0.488 0.62 0.68 0.244 0.224 0.316 0.884 0.912 ...
anova(modelo, test="F")

Portanto, ambas as variáveis são singificativas. Como ambas possuem apenas dois níveis, podemos parar por aqui e fazer o gráfico. Novamente, para fazer o gráfico, vou criar um dataframe com as informações de interesse, médias e erro padrão, além da identidade das variáveis. Aqui vou lançar mão de uma função muito interessante contida no pacote library(Rmisc), chamada summarySE. Nessa função preciso sempre informar quem são os dados, quem é a variável a ser medida (y) e quem são as variáveis agrupadoras (x).

library(Rmisc)
## Loading required package: plyr
(dados.plot <- summarySE(dados, measurevar = "Peso", groupvars = c("Sexo", "Turma")))

Vejam que tenho os dados agora, e posso prosseguir com o ggplot2. A explicação de como incorporar os dados será feita em sala.

ggplot(dados, aes(y=Peso, x=Turma, group=Sexo))+
  geom_jitter(aes(y=Peso, x=Turma, fill=Sexo), position = position_dodge(.9), pch=21)+
  geom_point(data= dados.plot, position = position_dodge(.9), size=3, pch=21, col="black", fill="grey50")+
  geom_errorbar(data= dados.plot, aes(ymax=Peso+se, ymin=Peso-se), position = position_dodge(.9), width=.3, col="grey10")

Pronto, finalizamos esse módulo com essas duas análises bastante utilizadas. Daqui por diante no curso, seguiremos compreendendo melhor o uso dos Modelos Lineares Generalizados.