Introdução ao R - Módulo IV
Autor: Pablo Fonseca
04 de Junho de 2018
Parte 1
Estatística descritiva utilizando o R
Nesse módulo iremos abordar algumas estratégias de análise de dados voltadas à análise de estatísticas descritivas de uma amostra. Logo em seguida, iremos abordar a criação de gráficos e ferramentas de visualização de dados.
Nesta aula iremos trabalhar alguns comandos e dicas para estatística descritiva usando o R.
A base de dados utilizada nesse módulo como exemplo será a base: “Modulo4_DataBase_CrohnD.txt” que está contida na pasta do DropBox do curso. Essa base de dados contém as informações de 117 pacientes afetados pela Doença de Crohn. Nessa base de dados temos o número de reações adversas (coluna nrAdvE) que cada paciente (coluna ID) apresentou em resposta ao tratamento de um novo fármaco. Existem três grupos nessa amostra, um grupo de individuos que possuem a doença de Crohn e receberam o fármaco número 1 (d1), um grupo de individuos que possuem a doença de Crohn e receberam o fármaco número 2 (d2) e um grupo de indivíduos que receberam placebo (placebo) (coluna treat). Além disso, existem algumas informações adicionais, como: índice de massa corporal (coluna BMI), altura em centímetros (coluna height), país de origem (coluna country), sexo (coluna sex), idade em anos (coluna age) e peso (coluna weight).
Mais informações sobre essa base de dados podem ser encontradas aqui: https://vincentarelbundock.github.io/Rdatasets/doc/robustbase/CrohnD.html
Em seguida iremos importar essa base de dados para o nosso ambiente de trabalho no R.
#Determinando o diretório de trabalho
setwd("E:/Doutorado/curso_R/Aula4/")
#Importando a base de dados
crohn<-read.table("Modulo4_DataBase_CrohnD.txt",h=T, sep="\t")
#Verificando a importação
head(crohn)## ID nrAdvE BMI height country sex age weight treat
## 1 19908 4 25.22 163 c1 F 47 67 placebo
## 2 19909 4 23.80 164 c1 F 53 64 d1
## 3 19910 1 23.05 164 c1 F 68 62 placebo
## 4 20908 1 25.71 165 c1 F 48 70 d2
## 5 20909 2 25.95 170 c1 F 67 75 placebo
## 6 20910 2 28.70 168 c1 F 54 81 d1#Verificando as dimensões da base de dados
dim(crohn)## [1] 117 9#Verificando a classe de cada uma das varáveis contidas na base de dados
str(crohn)## 'data.frame': 117 obs. of 9 variables:
## $ ID : int 19908 19909 19910 20908 20909 20910 21908 21909 21910 21911 ...
## $ nrAdvE : int 4 4 1 1 2 2 3 0 1 0 ...
## $ BMI : num 25.2 23.8 23.1 25.7 25.9 ...
## $ height : int 163 164 164 165 170 168 161 168 154 157 ...
## $ country: Factor w/ 2 levels "c1","c2": 1 1 1 1 1 1 1 1 1 1 ...
## $ sex : Factor w/ 2 levels "F","M": 1 1 1 1 1 1 1 1 1 1 ...
## $ age : int 47 53 68 48 67 54 53 53 47 58 ...
## $ weight : int 67 64 62 70 75 81 69 74 76 82 ...
## $ treat : Factor w/ 3 levels "d1","d2","placebo": 3 1 3 2 3 1 1 3 2 3 ...Note como a importação ocorreu como o esperado. A base de dados é composta por 117 pacientes e nove variáveis, já descritas anteriormente.
As análises descritivas que serão realizadas a seguir como exemplos serão baseadas nos valores referentes ao número de reações adversas que cada paciente apresentou durante o tratamento. Esses valores estão contidos na coluna “nrAdvE”.
Antes de iniciar a avaliação de medidas como média e mediana, iremos aplicar uma função muito interessante para a descrição do número de observações em diferentes classes dentro do nosso conjunto de dados.
A função table() contabiliza o número de observações contidos em cada uma das classes de variáveis do tipo fator ou de caracteres. A utilização da função é feita de maneira bem simples, como argumentos, devemos informar somente a(as) coluna(s) que desejamos contabilizar o número de observações.
A seguir, vamos contabilizar o número de pacientes de cada sexo, de cada país e que foi submetido a cada tipo de tratamento.
#Verificando o número de pacientes de cada sexo
table(crohn$sex)##
## F M
## 100 17#Verificando o número de pacientes de cada país
table(crohn$country)##
## c1 c2
## 78 39#Verificando o número de pacientes submetidos a cada tipo de tratamento
table(crohn$treat)##
## d1 d2 placebo
## 39 39 39Note que o output é composto por uma estrutura bem simples. Nele podemos notar a classe da variável e logo abaixo o número de observações. Uma combinação muito útil de funções envolvendo a função table() é combiná-la com a função as.data.frame() dessa maneira, nosso output será transformado em uma versão tabular. Vamos aplicar essa combinação utilizando como exemplo as classes de tratamento.
#Criando a versão tabular do comando table
table_ntreat<-as.data.frame(table(crohn$treat))
#Verificando se o output trata-se de uma Data frame
is.data.frame(table_ntreat)## [1] TRUE#Verificando o conteúdo da versão tabular
table_ntreat## Var1 Freq
## 1 d1 39
## 2 d2 39
## 3 placebo 39Note que o output realmente trata-se de uma versão tabular, ou seja, uma Data frame. Além disso, note como o output na versão tabular é constituído. Na primeira coluna (V1), estão contidas as classes das variáveis (d1, d2 e placebo). Na segunda coluna (Freq) temos o número de observações. Note que esses nomes são os nomes que a função determina para o output. Contudo, eles podems er alterado utilizando a função colnames().
Acima foi abordado a utilização da função table() para uma única coluna.Contudo, muitas vezes desejamos identificar o número de observações em função de mais de um grupo. Por exemplo, o número de participantes de cada sexo que recebeu cada um dos tratamentos. Para obter esses valores podemos executar o seguinte comando:
#Criando a versão tabular do output da função table
table.sexo.treat<-as.data.frame(table(crohn$sex,crohn$treat))
#Verificando o conteúdo do output
table.sexo.treat## Var1 Var2 Freq
## 1 F d1 34
## 2 M d1 5
## 3 F d2 34
## 4 M d2 5
## 5 F placebo 32
## 6 M placebo 7Note que como argumento da função table() foram informadas as colunas “sex” e “treat”. Além disso, agora, como colunas da Data frame gerada temos “V1”, “V2” e “Freq”. Com isso, conseguimos verificar o número de obervações por cada combinação de variáveis. Por exemplo, na base de dados existem somente 5 pacientes do sexo masculino tratados com o fármaco “d1”.
Uma alternativa muito útil da função table() é a a função prop.table(). com essa função podemos obter a proporção de cada uma das classes presentes na nossa base de dados. Vamos utilizar o mesmo exemplo acima para calcularmos os proporções de cada combinação entre o número de participantes de cada sexo que recebeu cada um dos tratamentos. Para isso, o output da função table() deve ser informado como argumento da função prop.table().
#Obtendo as proporções entre as combinações das classes
prop.table(table(crohn$sex,crohn$treat))##
## d1 d2 placebo
## F 0.29059829 0.29059829 0.27350427
## M 0.04273504 0.04273504 0.05982906Observe que de maneira geral, a proporção de pacientes do sexo masculino em todos os tratementos é menor. Além disso, note que os números presentes no output possuem uma longa sequência nas casas decimais. Em algusn momentos, apesar de não interferir em nada na representatividade do número, desejamos valores com um número de casas decimais reduzidos. Para isso, podemos utilizar a função round(). Com essa função podemos arredondar o número de casas decimais de qualquer número.
#Obtendo as proporções entre as combinações das classes
round(prop.table(table(crohn$sex,crohn$treat)),2)##
## d1 d2 placebo
## F 0.29 0.29 0.27
## M 0.04 0.04 0.06Note que no comando acima informamos como argumentos os valores criados à partir da função prop.table() e o valor 2. Com isso, informamos à função que desejamos arredondar todos os valores gerados pela função prop.table() em duas casas decimais.
Tente executar os comandos acima utilizando combinações entre diferentes classes de variáveis.
Um conjunto de funções muito úteis para estatísticas descritivas é composto pelas funções: min(), max(), mean(), median(), sd(). Essas funções apresentam os valores mínimo, máximo, a média, a mediana e o desvio-padrão de um conjunto de valores, respectivamente.
A seguir, vamos aplicar essas funções aos valores contidos na coluna “nrAdvE”.
#Obtendo o valor mínimo
min(crohn$nrAdvE)## [1] 0#Obtendo o valor máximo
max(crohn$nrAdvE)## [1] 12#Obtendo a média
mean(crohn$nrAdvE)## [1] 2.034188#Obtendo a mediana
median(crohn$nrAdvE)## [1] 1#Obtendo o desvio padrão
sd(crohn$nrAdvE)## [1] 2.728937Note que para todos os argumentos acima foi informado como argumento somente o conjunto de valores a ser analisado. Contudo, existe um detalhe importante para todas essas funções. Caso no conjunto de dados exista um missing value (NA), o resultado dessa função será NA. Vamos criar um objeto contendo 10 números aleatórios contidos no intervalo de 0 a 20 e um NA na décima primeira posição.
#Criando o objeto com elementos aleatório. Contudo, vamos utilizar a seed 123 para obtermos os mesmos valores
set.seed(123)
objA<-c(sample(0:20,10),NA)
#Verificando o conteúdo do objeto
objA## [1] 6 15 7 19 17 0 18 12 14 5 NA#Calculando a média desse objeto
mean(objA)## [1] NANote que o output da função mean() foi NA. Isso é devido à presença do NA na décima primeira posição do objeto. Nós poderíamos criar um novo objeto removendo os missing-values contido. Contudo, isso demanda tempo e esforços que não são necessários. Para isso, todo o conjunto de funções descrito acima aceita como argumento o comando na.rm=. Esse comando informa à função para que os NA sejam ignorados no cálculo.
#Criando o objeto com elementos aleatório. Contudo, vamos utilizar a seed 123 para obtermos os mesmos valores
set.seed(123)
objA<-c(sample(0:20,10),NA)
#Verificando o conteúdo do objeto
objA## [1] 6 15 7 19 17 0 18 12 14 5 NA#Calculando a média desse objeto
mean(objA, na.rm=T)## [1] 11.3Teste as demais funções utilizando o mesmo conjunto de dados e o argumento na.rm=.
Uma vez apresentado o conjunto de funções acima para obtermos alguns valores muito úteis para estatística descritiva, iremos aprender um comando que agrupa quase todas essas funções em uma única função. Essa função é a summary(). A função summary() apresenta como output o valor máximo, mínimo, a média, a mediana, o primeiro quartil e o terceito quartil. Note que somente o desvio-padrão não é obtido diretamente por essa função.
A seguir vamos aplicar a função summary() para o conjunto de dados contidos na coluna “nrAdvE” da Data frame “crohn”.
#Executando a função summary
summary(crohn$nrAdvE)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 1.000 2.034 3.000 12.000Note que o output da função segue a seguinte ordem: Min., 1st Qu., Median, Mean, 3rd Qu. e Max..
A função summary() consegue lidar bem com a presença de missing-values. Nesse caso, os valores apresentados no output da função já irão desconsiderar o NA. Porém, uma coluna estará presente no output contendo o número de missing-values contidos nos conjuntos (NA’s). Vamos utilizar o mesmo conjunto de dados criados anteriormente no objA.
#Executando a função summary
summary(objA)## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 0.00 6.25 13.00 11.30 16.50 19.00 1Nos comandos acima aprendemos alguns comandos básicos de estatística descritiva para a avaliação de um conjunto de dados sem separação por grupos. Porém, em diversos casos temos interesse de obter essas mesmas estatísticas para diferentes grupos de indivíduos em um mesmo conjunto de dados. Para isso, podemos utilizar uma função muito interessante, a função aggregate(). Os argumentos básicos da função aggregate() são: x=, by= e FUN. Esses três argumentos recebem a Data frame contendo o conjunto de dados a ser analisado, o conjunto de dados que agrupa as observações em diferentes grupos e a função a ser aplicada a esse conjunto de dados, respectivamente.
Por exemplo, vamos calcular a média de eventos adversos entre homens e mulheres:
#Obtendo a média de eventos adversos entre homens e mulheres
out_mean<-aggregate(x=crohn[,2], by=list(crohn$sex), mean)
#Avaliando o resultado da função aggregate
out_mean## Group.1 x
## 1 F 1.920000
## 2 M 2.705882#Verificando qual é a estrutura desse output
str(out_mean)## 'data.frame': 2 obs. of 2 variables:
## $ Group.1: Factor w/ 2 levels "F","M": 1 2
## $ x : num 1.92 2.71DETALHE: Note que o argumento by= deve ser informado como uma lista. Além disso, o argumento “crohn[,2]” informa que o conjunto de dados a ser analisado está inserido na segunda coluna da Data frame, no caso, a coluna “nrAdvE”.
Note que o output da função acima é uma Data frame contendo duas colunas. A primeira coluna (“Group.1”) contém cada um dos grupos (F e M) e a segunda coluna (“x”) contém os valores das médias para cada grupo.
Uma outra forma innteressante de se utilizar a função aggregate() é calculando os valores para a respectiva função escolhida para mais de um conjunto de valores por vez. Por exemplo, vamos calcular a média do número das reações adversar e da idade para cada sexo.
#Obtendo a média de eventos adversos e idade entre homens e mulheres
out_mean2<-aggregate(x=crohn[,c(2,7)], by=list(crohn$sex), mean)
#Avaliando o resultado da função aggregate
out_mean2## Group.1 nrAdvE age
## 1 F 1.920000 54.72000
## 2 M 2.705882 54.29412Com esse novo comando, onde as colunas 2 e 7 foram informadas no argumento contendo os dados a serem analisado, o output é composto por 3 colunas: A coluna contendo as classes avaliadas (F e M), a média para o número de reações adversas (nrAdvE) e a média da idade em cada grupo (Age).
Assim como podemos avaliar mais de uma variável por vez com a função aggregate(), podemos avaliar mais de um grupo. Como mostramos no uso da função table(), em alguns casos podemos ter o interesse de avaliar os valores da estatística descritiva em grupos compostos pela combinação de duas ou mais variáveis classificatórias em nosso conjunto de dados. Por exemplo, vamos calcular a mediana de reações adversos para homens e mulheres que foram submetidos a cada tipo de tratamento. Para isso, também utilizaremos a função aggregate().
#Calculando a mediana do número de reações adversas para homens e mulheres em função do tipo de tratamento
mean.by.groups<-aggregate(x=crohn[,"nrAdvE"], by=list(crohn$sex, crohn$treat), mean)
#Observando os valores obtidos
mean.by.groups## Group.1 Group.2 x
## 1 F d1 1.500000
## 2 M d1 2.000000
## 3 F d2 2.117647
## 4 M d2 2.200000
## 5 F placebo 2.156250
## 6 M placebo 3.571429Dentro da estatística descritiva existe uma medida que é amplamente utilizada e auxilia na avaliação da disperção dos dados em função da média de um determinado grupo. Essa medida é o z-score. O z-score pode ser facilmente calculado no R por meio da função scale().
O argumento básico a ser informado à função scale() é o conjunto de dados ao qual desejamos calcular os valores de z-score. No comando abaixo vamos calcular um z-score para o número de reações adversar observado para cada paciente, sem levar em considerção nenhum dos grupos. Em seguida, vamos inserir esse conjunto de dados em uma nova coluna de nossa base de dados. Essa coluna será nomeada como “z_nrAdvE”.
#Calculando os valores de z-score para cada medida contida em nrAdVE
crohn$z_nrAdvE<-scale(crohn$nrAdvE)
#Apresentando as seis primeiras linhas da Data frame
head(crohn)## ID nrAdvE BMI height country sex age weight treat z_nrAdvE
## 1 19908 4 25.22 163 c1 F 47 67 placebo 0.72035814
## 2 19909 4 23.80 164 c1 F 53 64 d1 0.72035814
## 3 19910 1 23.05 164 c1 F 68 62 placebo -0.37897102
## 4 20908 1 25.71 165 c1 F 48 70 d2 -0.37897102
## 5 20909 2 25.95 170 c1 F 67 75 placebo -0.01252797
## 6 20910 2 28.70 168 c1 F 54 81 d1 -0.01252797Observe que na Data frame “crohn” a nova coluna, “z_nrAdvE”, contém os valores obtidos para o z-score do número de reações adversas de todos os pacientes contidos na tabela.
Uma abordagem frequentemente utilizada para o cálculo de z-score é o cálculo do mesmo em função de grupos presentes na amostra. Nesse caso, o desvio é calculado em função da média de cada grupo que cada mediada pertence. Na Data frame que estamos utilizando como exemplo podemos realizar esse procedimento para calcular o z-score no número de reações adversas de cada paciente em função dos grupos criados pela variável tratamento. Esse procedimento por ser feito de diversas formas. Aqui, de modo a apresentar mais formas de se trabalhar com diferentes grupos dentro de funções básicas, iremos executar esse cálculo por meio da função ave(). A função ave() recebe como argumentos o vetor contendo os valores a serem utilziados para o cálculo do z-score, a variável contendo o grupo em que cada paciente está inserido, e a função a ser utilizada, no caso, a função scale(). Ao término do processo, iremos criar uma nova coluna, denominada “z_nrAdvE_treat” para inserir os valores.
#Calculando os valores de z-score para cada medida contida em nrAdVE em função dos grupos presentes em treat
crohn$z_nrAdvE_treat<-ave(x=crohn[,"nrAdvE"], crohn$treat, FUN = scale)
#Apresentando as seis primeiras linhas da Data frame
head(crohn)## ID nrAdvE BMI height country sex age weight treat z_nrAdvE
## 1 19908 4 25.22 163 c1 F 47 67 placebo 0.72035814
## 2 19909 4 23.80 164 c1 F 53 64 d1 0.72035814
## 3 19910 1 23.05 164 c1 F 68 62 placebo -0.37897102
## 4 20908 1 25.71 165 c1 F 48 70 d2 -0.37897102
## 5 20909 2 25.95 170 c1 F 67 75 placebo -0.01252797
## 6 20910 2 28.70 168 c1 F 54 81 d1 -0.01252797
## z_nrAdvE_treat
## 1 0.5044258
## 2 1.0472493
## 3 -0.4474745
## 4 -0.4261340
## 5 -0.1301744
## 6 0.1874025Compare os valores contidos nas duas colunas geradas por meio do uso da função scale(), note que os valores diferem-se substancialmente. Esse é um procedimento muito útil para a padronização de dados para análises estatísticas futuras.
Uma outra forma de se utilizar o z-score, que é frequentemente aplicada, é o cálculo do mesmo dentro de faixas de valores. Por exemplo, dentro de faixas de idade. Iremos desenvolver um exemplo para calcularmos o z-score no número de reações adversas de cada paciente em função da média de três faixas etárias. Para isso, iremos criar três novos grupos. O primeiro grupo irá ser constituído pelos pacientes com idades entre 40 até 50 anos. O segundo grupo por paciente com idades entre 50 até 60 anos. O terceiro grupo, por pacientes com idade superior a 60 anos. Essa nova classificação de grupos será inserida em uma nova coluna, que será denominada “age_interval”.
#Preenchendo a nova coluna com missing values que serão futuramente substituídos
crohn$age_interval<-rep(NA,nrow(crohn))
#Criando os grupos cbaseado nas faixas etárias
crohn[which(crohn$age>=40 & crohn$age<50),"age_interval"]<-"40-50"
crohn[which(crohn$age>=50 & crohn$age<60),"age_interval"]<-"50-60"
crohn[which(crohn$age>60),"age_interval"]<-">60"
#Calculando o z-score em função da faixa etária
crohn$z_nrAdvE_age_int<-ave(x=crohn[,"nrAdvE"], crohn$age_interval, FUN = scale)
#Avaliando o resultado
head(crohn)## ID nrAdvE BMI height country sex age weight treat z_nrAdvE
## 1 19908 4 25.22 163 c1 F 47 67 placebo 0.72035814
## 2 19909 4 23.80 164 c1 F 53 64 d1 0.72035814
## 3 19910 1 23.05 164 c1 F 68 62 placebo -0.37897102
## 4 20908 1 25.71 165 c1 F 48 70 d2 -0.37897102
## 5 20909 2 25.95 170 c1 F 67 75 placebo -0.01252797
## 6 20910 2 28.70 168 c1 F 54 81 d1 -0.01252797
## z_nrAdvE_treat age_interval z_nrAdvE_age_int
## 1 0.5044258 40-50 0.5641145
## 2 1.0472493 50-60 1.2431607
## 3 -0.4474745 >60 -0.5037838
## 4 -0.4261340 40-50 -0.4615482
## 5 -0.1301744 >60 -0.2032812
## 6 0.1874025 50-60 0.2685227Note como o conjunto de comandos acima foi capaz de criar grupos em função da faixa etária de maneira eficaz. Novamente, compare os valores de z-score e veja como eles diferem-se entre si.
Testes de associação utilizando o R
O foco desse módulo não será sobre testes de significância ou comparação de média entre grupos. Contudo, iremos falar brevemente sobre os teste de qui-quadrado, teste exato de fisher, teste-t e ANOVA utilizando o R.
É importante ressaltar que o R possui praticamente toda a gama de testes estatísticos possível. A busca por informações sobre esses testes e a forma de implementação no R pode ser fácilmente realizada por meio de fóruns de discussão ou nos manuais dos pacotes que implementaram esses testes. Dessa forma, é importante relembrar o que foi dito na primeira aula: Busque os fóruns de discussão online e não deixe de utilizar o google para encontrar o funcionamento básico das funções de interesse.
Qui-quadrado e Teste exato de Fisher
O qui-quadrado e o teste exato de Fisher são testes de associação que comparam a frequência de observações de um determinado evento entre diferentes grupos. No R, o qui-quadrado e o teste exato de Fisher são calculados utilizando as funções chisq.test() e fisher.test(), respectivamente. Para essas funções, devemos passar como argumentos os conjuntos de dados, que agrupam nossas observações em diferentes grupos, que temos interesse de comparar. Utilizando a Data frame “crohn” como exemplo, iremos calcular se existem diferenças siginificativas entre a frequência de homens e mulheres dentro de cada grupo de país avaliado. Abaixo iremos apresentar uma nova função, a função addmargins(). Essa função utiliza a tabela de contigência criada pela função table() e cria o somatório de cada grupo nas margens. Essa é uma função interessante para a verificação dos dados. Detalhe, nós não iremos inserir as margens na tabela de contigência que será utilziada para o teste do qui-quadrado.
#Obtendo a tabela de contigência
table.country.treat<-table(crohn$country,crohn$sex)
#Verificando a tabela criada
table.country.treat##
## F M
## c1 65 13
## c2 35 4#Veriricando o sometório das margens da tabela
addmargins(table.country.treat)##
## F M Sum
## c1 65 13 78
## c2 35 4 39
## Sum 100 17 117#Calculando o qui-quadrado
chisq.test(crohn$country, crohn$sex)##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: crohn$country and crohn$sex
## X-squared = 0.42154, df = 1, p-value = 0.5162#Calculando o teste exato de fisher
fisher.test(crohn$country, crohn$sex)##
## Fisher's Exact Test for Count Data
##
## data: crohn$country and crohn$sex
## p-value = 0.4169
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.1267055 2.0455382
## sample estimates:
## odds ratio
## 0.5739733Note o output das funções chisq.test() e fisher.test(). O output da função chisq.test() informa que o teste de qui-quadrado foi realizado aplicando a correção de continuidade de Yates, utilizada quando existem caselas com valores iguais ou inferiores a 5 observações. O uso dessa correção é padrão para essa função. Caso seja de interesse não aplicar tal correção, devemos utilizar o argumento “correct=FALSE” dentro da função chisq.test(). Além disso, note que são apresentados os dados avaliados logo após “data:” e os valores de qui-quadrado (X-squared), os valores de graus de liberdado (df) e o valor de p (p-value) do teste. Note que não há diferenças siginificativas entre as observações. O teste exato de Fisher é uma boa alternativa para quando possuímos um conjunto de dados onde existem caselas com um baixo número de observações. Repare no output do teste exato de Fisher. O output apresenta os dados utilizados, o valor de p, uma descrição do teste de hipotese utilizado (no caso, h0 é definido como uma razão de change igual a 1), os limites inferiores e superiores para um intervalo de confiança de 95% dos dados e a razão de chance calculada entre os grupos.
DICA: A distribuição utilizada para o cálculo do p-value na função chisq.test() pode ser simulada levando-se em consideração os dados analisados. Isso pode ser feito por meio da utilização dos argumentos “simulate.p.value=” e “B=”, onde, o primeiro recebe uma resposta lógica (TRUE) no caso da escolha pela simulação e o segundo recebe o número de replicatas para as simulações realizadas (o padrão pré-definido para B é 2000). Para mais informações utilize o comando ?chisq.test().
Teste-t
O teste-t é comumente utilizado para comparar a média entre grupos. No R, o teste-t é executado pela função t.test(). Um detalhe importante do uso da função t.teste(), na verdade, da utilização do teste-t, é que precisamos saber se a variância entre os grupos é igual. Para isso, utilizamos o teste-F, implementado na função var.test(), no R. Para exemplificarmos o funcionamento dessas funções, iremos comparar o número de reações adversas entre cada tratamento da Data frame “crohn”.
#Verificando se a variância é igual entre os grupos
var.test(crohn$nrAdvE ~ crohn$treat)## Error in var.test.formula(crohn$nrAdvE ~ crohn$treat): grouping factor must have exactly 2 levelsDETALHE: O sinal “~” é constantemente utilizado no R para representar a dependência entre duas variáveis. Assim, no comando acima estamos representando que iremos avaliar a igualdade de variância de “nrAdvE” em função dos grupos contidos em “treat”.
Note que o output da função var.test() informa que o número de grupos deve ser exatamente igual a 2. Isso é uma mensagem de erro referente a uma limitação do teste-F. Só podemos comparar a variância de dois grupos. Por exemplo, se tentarmos realizar o mesmo teste, porém, com os grupos baseados no sexo dos paciente, não encontraremos essa mensagem de erro.
#Verificando se a variância é igual entre os grupos
var.test(crohn$nrAdvE ~ crohn$sex)##
## F test to compare two variances
##
## data: crohn$nrAdvE by crohn$sex
## F = 0.67306, num df = 99, denom df = 16, p-value = 0.2404
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 0.2808052 1.3061147
## sample estimates:
## ratio of variances
## 0.673063Note que o output é semelhante ao observado no teste exato de Fisher. A H0 desse teste é definida como a razão entre as variâncias dos grupos ser igual a 1. Como observamos um p-value maior que 0.05 (utilziando 0.05 como o nosso alfa, o limiar de significância), não rejeitamos h0. Ou seja, as variâncias para o número de reações adversas entre os sexos são iguais.
Existem algumas alternativas para a comparação de variâncias entre mais de 2 grupo (no caso, nosso grupo tratamento possui 3 grupos). Uma das mais utilziadas é o teste de Bartlett. No R, esse teste está implementado na função bartlett.test(). A seguir vamos comparar as variâncias entre o número de reações adversas para cada um dos tratamentos.
DETALHE: Iremos considerar que a distribuição do dado é normal. Iremos falar um pouco sobre testes de normalidade na próxima aula. No caso do dado não possuir uma distribuição normal, o Teste de Levene é uma alternativa mais robusta do que o Teste de Barlett.
#Verificando se a variância é igual entre os grupos
bartlett.test(crohn$nrAdvE ~ crohn$treat)##
## Bartlett test of homogeneity of variances
##
## data: crohn$nrAdvE by crohn$treat
## Bartlett's K-squared = 3.5146, df = 2, p-value = 0.1725Note que o valor de p obtido é maior que 0.05. Portanto, consideramos que as variâncias entre os grupos são iguais. Com isso, podemos definir de maneira adequada os argumentos para o teste-t. Novamente, o teste-t consegue lidar com comparações par a par. Portanto, só podemos comparar duas amostras por vez.
#Aplicando o teste-t para o número de reações adversas ente os tratamentos, utilizando igualdade de variância entre os grupos
#Droga 1 x placebo
t.test(crohn[which(crohn$treat=="d1"),"nrAdvE"], crohn[which(crohn$treat=="placebo"),"nrAdvE"], var.equal=T)##
## Two Sample t-test
##
## data: crohn[which(crohn$treat == "d1"), "nrAdvE"] and crohn[which(crohn$treat == "placebo"), "nrAdvE"]
## t = -1.3491, df = 76, p-value = 0.1813
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.0953702 0.4030626
## sample estimates:
## mean of x mean of y
## 1.564103 2.410256#Droga 2 x placebo
t.test(crohn[which(crohn$treat=="d2"),"nrAdvE"], crohn[which(crohn$treat=="placebo"),"nrAdvE"], var.equal=T)##
## Two Sample t-test
##
## data: crohn[which(crohn$treat == "d2"), "nrAdvE"] and crohn[which(crohn$treat == "placebo"), "nrAdvE"]
## t = -0.42794, df = 76, p-value = 0.6699
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -1.594758 1.030655
## sample estimates:
## mean of x mean of y
## 2.128205 2.410256#Droga 1 x Droga 2
t.test(crohn[which(crohn$treat=="d1"),"nrAdvE"], crohn[which(crohn$treat=="d2"),"nrAdvE"], var.equal=T)##
## Two Sample t-test
##
## data: crohn[which(crohn$treat == "d1"), "nrAdvE"] and crohn[which(crohn$treat == "d2"), "nrAdvE"]
## t = -0.99962, df = 76, p-value = 0.3207
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -1.6880377 0.5598326
## sample estimates:
## mean of x mean of y
## 1.564103 2.128205Note que como argumento para a função t-test() devemos informar os conjuntos de observações para cada grupo de maneira separa. Por isso tivermos que realizar o processo de filtragem dos dados utilziando o sistema de busca e indexação já trabalhado ensse curso. Além disso, note que o argumento “var.equal=” teve como valor informado TRUE (T), isso indica que a variância entre os grupos é igual.
Agora, note o output do teste-t. Logo após “data:” é apresentado o conjunto de dados avaliado, representado pelos dois primeiros argumentos informados à função t.teste(). Em seguida, temos os valores da distribuição t, os graus de liberdade e os valores de p para o teste realizado. Em seguida é apresentado o teste de hipotese realizado, nesse caso, H0 é definido como a igualdade entre as médias representado pela diferença entre as médias de cada grupo ser igual a zero. Em seguida, são apresentados os limites inferiores e superiores para o intervalo de confiança de 95% dos dados. Por fim, a média de cada grupo é apresentada.
Note que nos três testes, não houve diferença siginificativa entre as médias dos grupos nos três testes (P>0.05).
ANOVA
Uma alternativa rápida e simples para a comparação da média de mais de dois grupos é o uso da ANOVA. a ANOVA está implementada no R na função aov().
#Calculando a ANOVA do número de reações adversas em função dos tipos de tratamento
anova_out<-aov(crohn$nrAdvE ~ crohn$treat)
#Obtendo os valores do teste estatístico da comparação das médias entre os grupos
summary(anova_out)## Df Sum Sq Mean Sq F value Pr(>F)
## crohn$treat 2 14.5 7.239 0.972 0.382
## Residuals 114 849.4 7.451IMPORTANTE: Note que devemos utilizar a função summary() no output da função oav() de modo a obter os valores de p para a comparação entre as médias dos grupos.
Observe que, assim como no teste-t, o valor de p obtido para a comparação entre as médias dos grupos não é significativo (P>0.05). É importante ressaltar que o valor de p está contido na coluna “Pr(>F)”, que é a probabilidade de se encontrar um valor superior ao valor da estatística F encontrada para o teste realizado.
É importante ressaltar que a ANOVA, da forma que foi calculada, compara se existe diferença entre as médias de pelo menos um dos grupos. Portanto, caso encontrassemos um valor de p inferior a 0.05, não seríamos capazes de definir entre quais grupos observamos essa diferença em um primeiro momento. Para isso, seria necessário executarmos um post-hoc. Nessa aula não iremos abordar o uso de post-hoc. Contudo, um dos testes mais utilzados para post-hoc é o o teste de Tukey. Mais informações sobre esse teste e sua implementação no R podem ser encontrados aqui: https://www.r-bloggers.com/anova-and-tukeys-test-on-r/
Nessa aula abordamos uma série de funções muito úteis para a execução de estatística descritiva utilizando o R. Na próxima aula iremos abordar a utilização de funções que auxiliam na visualização gráfica dos nossos dados em forma de gráficos (gráficos de pontos, histogramas, box-plot, gráficos de barra, etc.). A utilização simultânea das funções apresentadas nessa aula com as funções de visualização gráfica que serão apresentadas na próxima aula tornam o processo de obtenção da estatística descritiva dos dados muito mais completo e confiável.
Parte 2
Gráficos utilizando o R
Uma das grandes vantagens de se utilizar o R é a capacidade que o mesmo possui de produzir gráficos e figuras de maneira extremamente eficiente e representativa. Nesse módulo, iremos apresentar os conceitos básicos que envolvem a criação e edição de gráficos no R. além disso, iremos abordar alguns tipos específicos de gráficos que são constantemente utilizados em rotinas de pesquisa, como por exemplo, Gráfico de barras, Box-plots, Histogramas, Gráfico de pizza e Heatmap.
Nessa aula, iremos utilizar a base de dados apresentada na aula anterior para criarmos a maior parte dos exemplos mostrados. A base de dados é nomeada “Modulo4_DataBase_CrohnD.txt” e está disponível na pasta do Dropbox do curso.
Contudo, para iniciarmos, iremos trabalhar com um conjunto de dados mais simples. Com isso, poderemos focar nos argumentos principais utilizados na função mais básica para a criação de gráficos (plots) utilizando o R: a função plot().
Iremos criar dois objetos para introduzirmos a utilização da função plot(). O primeiro objeto será composto por uma sequência de números que inicia em 1 e vai até 10 sendo incrementado de 1 em 1. O segundo objeto será composto por uma sequência de números que inicia em 10 e vai até 100 sendo incrementado de 10 em 10.
#Criando os objetos contendo as sequências de números
objeto1<-seq(from=1, to=10, by=1)
objeto1## [1] 1 2 3 4 5 6 7 8 9 10objeto2<-seq(from=10, to=100, by=10)
objeto2## [1] 10 20 30 40 50 60 70 80 90 100O primeiro gráfico a ser criado será o mais simples possível, utilizando os argumentos essenciais para o funcionamento da função plot(). Os dois argumentos básicos da função plot são os argumentos “x=” e “y=”. Esses argumentos recebem como input os elementos a serem plotados no eixo-x e eixo-y, respectivamente. No exemplo a seguir, iremos plotar o objeto1 no eixo-x e o objeto2 no eixo-y.
#Plotando objeto1 e objeto2
plot(x=objeto1, y=objeto2)
Note a estrutura básica do gráfico criado. Ele é composto por pontos, os eixos x e y possuem uma sequência pré-determinada para a apresentação dos valores e o nome dos eixos é referente ao nome dos objetos utilizados nos argumentos.
Esse é um gráfico muito simples, apesar de apresentar todos os componentes necessários. Dessa forma, podemos aplicar uma série de melhorias no gráfico. A seguir, iremos trabalhar de uma maneira sequencial, introduzindo novos argumentos que irão realizar edições nos gráficos baseados em 3 grupos: edições envolvendo os eixos, edições envolvendo os caracteres plotados e edições envolvendo a estrutura do plot.
- Edições envolvendo os eixos
Partindo do plot criado anteriormente, iremos apresentar uma série de argumentos relacionados com a edição dos eixos presentes no gráfico. Com esses comandos iremos alterar algumas características do plot como o nome presente nos eixos, o tamanho da fonte utilizada nos eixos, inserir títulos nos gráficos, alterar o tamanho dos números apresentados nos eixos, o intervalos plotados nos eixos, etc.
Abaixo iremos alterar o nome dos eixos x e y e iremos inserir um título para o gráfico, por meio dos argumentos “xlab=”, “ylab=” e “main=”, respectivamente. Além disso, iremos aumentar o tamanho da fonte utilizada no nome dos eixos do plot por meio do argumento “cex.lab=”. Por fim, iremos aumentar a fonte dos números apresentados nos eixos do gráfico por meio do argumento “cex.axis=”.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2)
Note o gráfico gerado pelo comando acima e compare com o primeiro gráfico criado. O nome dos eixos foram alterados em função do nome informado para os argumentos “xlab=” e “ylab=”. Além disso, repare que o tamanho da fonte do título dos eixos e dos números contidos nos eixos foi aumentado. Esse aumento foi relizado por meio dos argumentos “cex.lab=” e “cex.axis=”. É importante ressaltar que os valores informados para esses argumentos funcionam de uma maneira referencial. O valor referencial para a plotagem das fontes é 1. Assim, valores maiores que 1 resultam em um aumento proporcional da fonte e valores menores que 1, resultam em uma diminuição proporcional da fonte. Por fim, note que um título foi inserido no gráfico, referente ao nome informado à função “main=”.
DICA: Procure alterar o nome dos eixos e o tamanho das fontes com valores à sua escolha.
Além das alterações realizadas acima, ainda podemos fazer algumas modificações referentes aos eixos. A seguir, vamos inserir mais alguns argumentos no comando utilizado previamente para gerar o gráfico. Nesse novo comando, além das modificações anteriores, iremos alterar os valores minímos e máximos presentes nas escalas dos eixos x e y por meio dos argumentos “xlim=” e “ylim=”, respectivamente.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,20), ylim=c(0,150))
Note que no gráfico apresentado os eixos tiveram seus intervalos de valores alterados. Agora, o eixo-y apresenta valores de 0 até 150 e o eixo-x apresenta valores de 0 a 20. Contudo, é importante notar que os valores presentes nos argumentos “x=” e “y=” não foram alterados. Portando a distribuição dos pontos dentro do gráfico não foi alterada. A utilização dos argumentos “xlim=” e “ylim=” deve ser feita de maneira bem precisa, de modo a não afetar as proporções dos dados (como aconteceu no gráfico gerado no comando anterior). Para os próximos gráficos iremos utilizar como valores mínimos e máximos para os eixos y e x c(0,100) e c(0,10), respectivamente.
- Edições envolvendo os caracteres plotados
A seguir iremos adicionar alguns argumentos que alteram as características dos caracteres plotados no gráfico.
O primeiro argumentos que iremos utilizar é o argumento “type=”. Esse argumento quando especificado dentro da função plot(), permite a seleção do tipo de caracter plotado. Por exemplo, caso seja necessário plotar os dados em linha ao invés de pontos, podemos informar a esse argumento o valor “l”. Pontos são os caracteres padrão da função plot, o valor que representa pontos para o argmento “type=” é “p”. No caso de ser necessário plotar pontos e linhas, podemos utilizar o valor “b”.
Mais informações sobre os valores aceitos pelo argumento “type=” podem ser encontrados utilizando o comando ?plot().
A seguir vamor transformar o nosso gráfico em um gráfico contando linhas e pontos.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b")
Note que pela simples adição do argumento “type=”b“”, mudamos a forma como os dados foram apresentados. Tente mudar a apresentação dos dados por meio de outros valores informados para o argumento “type=”.
Uma outra forma de alteração dos caracteres apresentados é por meio do argumento “pch=”. Esse argumento representa qual o caracter que será plotado no gráfico. Nos comandos anteriores fomos capazes de notar que o caracter padrão é um círculo sem preenchimento. Contudo, como apresentado na figura 1, vários caracteres podem ser utilizados na criação de gráficos no R. Esses caracteres recebem valores específicos que deverão ser informados para o argumento “pch=”.
Figura 1: Codificação para diferentes tipos de caracteres a serem plotados no R.Note que na figura as cores vermelho e amarelo foram escolhidas para apresentas os caracteres. Contudo, o padrão de apresentação no R é a cor preta.
Abaixo, iremos editar o gráfico anterior de maneira a apresentarmos um gráfico contendo linhas e triângulos preenchidos (pch=17) como caracteres representando os dados.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17)
Tente criar gráficos com diferentes tipos de pch.
Um outro conjunto de alterações envolvendo os caracteres presentes nos plots é a alteração do tamanho dos caracteres e sua cor. Essas alterações são feitas por meio dos argumentos “cex=” e “col=”, respectivamente. Abaixo iremos dobrar o tamanho dos caracteres (cex=2) e alterarmos a cor para vermelho (col=“red”).
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red")
Note como o tamanho dos critérios cresceu de maneira proporcional ao valor de referência (cex=1). Além disso, agora os caracteres são vermelhos.
Existe uma infinidade de cores possíveis de serem selecionadas no R. Para um melhor detalhamento das opções de cores no R, acesse o link a seguir para visualizar a palheta de cores disponível no R, assim como os nomes que devem ser informados ao argumento “col=” para obter essas cores: http://www.stat.columbia.edu/~tzheng/files/Rcolor.pdf
- Edições envolvendo a estrutura do plot
Esse terceiro grupo de alterações irá ser referente a alterações que podemos fazer na composição do plot. Seja na presença de margem em torno do plot, na sequência númerica contida nos eixos ou em parâmetros gráficos da janela do plot.
A primeira alteração que iremos abordar é a remoção da linha que compõem a margem em torno do plot. Essa é uma modificação puramente estética. Contudo, essa é uma das vantagens do R, você consegue editar praticamente todos os critérios estéticos do seus gráficos. Essa remoção das margens pode ser feito por meio do argumento “bty=. Para removermos as margens, devemos informar a esse argumento o valor”n“.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red", bty="n")
Em alguns casos, a eliminação das margens auxiliam na visualização dos dados. Contudo, como dito anteriormente trata-se mais deum critério estético e de preferencial pessoal.
Uma outro alteração na estrutura do dos plots é a alteração do padrão de incremento dos eixos x e y. Baseado na distribuição dos valores utilizados, o R seleciona uma escala de incremento para os eixos. No caso do exemplo que estamos utilizando, o eixo-y foi incrementado de 20 em 20 e o eixo-x de 2 em 2. Em alguns casos, para facilitar a interpretação do dados, podemos querer alterar esse padrão. Isso é possível por meio da utilização de uma outra função, complementar à função plot. Essa é a função axis(). Com essa função somos capazes de informar vários argumentos que irão guiar o padrão de preenchimento dos eixos do gráfico. Contudo, antes de criarmos novos eixos, precisamos desabilitar os eixos criados automaticamente pela função plot(). Essa desabilitação dos eixos é feito pelo argumento “axes=”. Para removermos os eixos padrão da função plot(), precisamos informar o valor lógico “FALSE” para o argumento “axes=”.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red", bty="n", axes=F)
Note que no ouput acima os eixos foram desabilitados. Dessa forma, somos capazes de criarmos novos eixos, seguindo o padrão desejados por nós, para o gráfico.
Uma vez retirados os eixos, a função axis() pode ser utilizada logo após a função plot() para criarmos os novos eixos. O funcionamento dessa função é bem simples. Os argumentos básicos a serem informados para a função axis() são: “side=”, “at=” e “labels=”. Esses argumentos determinam qual eixo será preenchido, qual a posição dos valores no eixo e a sequência de valores a ser preencida, respectivamente. O argumento “side=” recebe valores que vão de 1 a 4: 1= eixo inferior, 2= eixo lateral esquerdo, 3= eixo superior, 4= eixo lateral direito. Ou seja, partindo do eixo inferior (eixo-x), os valores aumentam até 4 seguindo o sentindo horário para os 4 lados do plot. Portanto, como desejamos determinar uma nova sequência para os eixos x e y, devemos executar a função axis() duas vezes, uma informando o valor 1 para o argumento *side=" e uma segunda informando o valor 2 para o argumento “side=”. Assim, iremos criar novas sequências para os eixos x e y, respectivamente.
Abaixo iremos incrementar o eixo x de 1 em 1 e o eixo y de 10 em 10.
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red", bty="n", axes=F)
#Criando o novo eixo-x contendo valores incrementados de 1 em 1
axis(side=1, at=seq(1,10,1), labels=seq(1,10,1), cex.axis=0.7)
#Criando o novo eixo-y contendo valores incrementados de 10 em 10
axis(side=2, at=seq(10,100,10), labels=seq(10,100,10), cex.axis=0.7)
Notem como agora os eixos estão sendo incrementados nos respectivos padrões selecionados. É importante ressaltar que eu optei por diminuir o tamanho da fonte dos valores nos eixos para que fosse possível notarmos todos os valores nos eixos. Contudo, podemos continuar mantendo uma fonte maior e resolvermos o problema da falta de espaço para que todos os valores apareçam no momento da exportação do gráfico.
A última alteração envolvendo a estrutura da janela do plot que será apresentada está relacionada com o número de plots presentes em uma única janela. No R, existem diversas formas de criar janelas contendo mais de um gráfico. Aqui, iremos iremos utilizar a função mfrow(). Essa função é na verdade um parâmetro gráfico que pode ser informada à função par(). A função par() é responsável por determinar uma série de parâmetros gráficos que serão aplicados a todos os plots criados após a sua utilização. Aqui, iremos focar somente na criação de janelas com mais de um plot. Para isso, devemos informar à função mfrow() o número de linhas e colunas que teremos em nossa janela, onde serão plotados os gráficos.
A seguir, iremos plotar o primeiro gráfico criado nessa aula ao lodo do último gráfico criado. Para isso, devemos informar que teremos uma linha e duas colunas em nossa janela.
#Determinando os parâmetros gráficos
par(mfrow=c(1,2))
#Plotando objeto1 e objeto2 como no primeiro gráfico
plot(x=objeto1, y=objeto2)
#Plotando o segundo gráfico
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red", bty="n", axes=F)
#Criando o novo eixo-x contendo valores incrementados de 1 em 1
axis(side=1, at=seq(1,10,1), labels=seq(1,10,1), cex.axis=0.7)
#Criando o novo eixo-y contendo valores incrementados de 10 em 10
axis(side=2, at=seq(10,100,10), labels=seq(10,100,10), cex.axis=0.7)
#Resetando os prâmetros gráficos anteriores
dev.off()## null device
## 1IMPORTANTE: Note que após determinaros novos parâmetros gráficos e executarmos todos os comandos envolvidos com essas alterações, o comando dev.off() foi utilizado. Esse comando precisa ser executado para que os parâmetros gráficos sejam “resetados” para os valores padrão.
Exportando os gráficos criados para uma pasta no computador
A seguir, iremos apresentar 3 formas de exportar um gráfico criado no ambiente para o seu computador. Para os exemplos a seguir, utilizaremos o último plot criado. Contudo, ao invés de utilizarmos cex.axis=0.7, iremos utilizar cex.axis=2. Nesse caso, iremos apresentar como fazer para que todos os números dos eixos tornem-se visíveis no output.
- Exportando gráficos utilizando a janela Plots
Na figura 2 são apresentados alguns dos passos para exportar um gráfico criado, por meio da aba Plot. Uma vez que o gráfico foi criado, ele será apresentado na aba plot, na janela inferior direita do R studio. Nessa janela, podemos notar a presença do botão “export”. Clicando nesse botão, nos é apresentado as opções de exportação do gráfico. As opções apresentadas são: “Save as Image…”, “Save as pdf…” e “Copy to clipboard”. A opção “Copy to clipboard” simplesmente realiza um “ctrl+c” da imagem, possibilitando colá-la em outro software, como o Power Point. Aqui, iremos focar nas opções apresentadas para salvar a figura como imagem.
Após a seleção de “Save as Image…”, uma nova janela é aberta. Nessa janela, no canto superior esquerdo, temos a opção de salvar como diferentes formatos de imagem. Nesse exemplo, iremos selecionar o formato “.png”. Após a seleção do formato, logo abaixo existe o botão “Directory”. Clicando nesse botão somos capazes de selecionar em qual diretório desejamos salvar o gráfico. O detalhe é que o diretório de trabalho atual já estará selecionado. Logo abaixo do botão “Directory”, existe a opção “File name”. Nessa opção devemos informar qual o nome que desejamos dar ao arquivo “.png” que será criado.
No canto superior direito estão localizadas as opções que irão nos permitir salvar o gráfico de uma maneira onde todos os números presentes nos eixos sejam visíveis, mesmo aumentado o tamanho da fonte. Essas opção são: “width” e “height”. Essas opções especificam a largura e altura do gráfico a ser salvo, especificamente. No exemplo contido na figura 2, determinamos que o gráfico possuirá 1200 por 1200 pixels. Por fim, clique em “save”, no canto inferior direito para salvar o gráfico com as opções selecionadas. Note que no output presente na pasta do Dropbox, “plot_modulo4_aula2_utilizando_janela.png”, todos números são visíveis.
Figura 2: Passo a passo para exportação do gráfico criado por meio da aba Plot.
- Exportando uma imagem no formato .png diretamente pelo terminal
Além da utilização da aba Plot para exportar um arquivo no formato “.png”, existe a posibilidade de realizar o mesmo procedimento por meio do terminal. Isso é possível por meio da função “png()”. Essa função deve ser executada anteriormente aos comandos que irão gerar o gráfico. A função “png()” recebe como argumentos básicos “file=”, “height=” e “width=”. O argumento “file=” recebe o nome com que o arquivo será salvo. É importante ressaltar que no nome informado ao argumento “file=”, devemos também adicionar a extensão “.png” ao nome. Os argumentos “width=” e “height=” recebem os valores referentes à largura e altura do arquivo. Além disso, existem outros dois argumentos que podem ser informados: “units=” e “res=”. Esses argumentos alteram a unidade dos valores informados para os argumentos “width=” e “height=” e a resolução da imagem criada, respectivamente. No exemplo abaixo não iremos abordar esses dois argumentos. Contudo, mais informações podem ser encontradas aqui: https://stat.ethz.ch/R-manual/R-devel/library/grDevices/html/png.html
#Criando o arquivo .png que irá receber o gráfico
png(file="plot_modulo4_aula2_utilizando_terminal.png", width=1200, height=1200)
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red", bty="n", axes=F)
#Criando o novo eixo-x contendo valores incrementados de 1 em 1
axis(side=1, at=seq(1,10,1), labels=seq(1,10,1), cex.axis=0.7)
#Criando o novo eixo-y contendo valores incrementados de 10 em 10
axis(side=2, at=seq(10,100,10), labels=seq(10,100,10), cex.axis=0.7)
#Resetando os prâmetros gráficos anteriores
dev.off()## png
## 2Note que assim como quando criamos uma janela contendo dois gráficos, a função dev.off() foi executada para “resetar” os parâmetros gráficos. Isso é de extrema importância, uma vez que o gráfico só será realmente salvo no arquivo “.png” se essa função for executada. Repare no arquivo “plot_modulo4_aula2_utilizando_terminal.png” que foi gerado com o comando acima e está presente na pasta do Dropbox. Esse arquivo é exatamente igual ao arquivo criado pela exportação por meio da aba Plot.
VANTAGEM: A vantagem de utilizar a função png(), ou outras semelhantes (bmp(), tiff(), jpeg(), etc.), é que essas funções permitem uma automatização do processo, não dependendo que tenhamos que executar o comando manualmente, assim interrompendo o andamento do script.
- Exportando arquivos no formato .pdf diretamente pelo terminal
Uma outra forma de exportarmos imagens e gráficos criados no ambiente e por meio de arquivos do tipo “.pdf”. Esse procedimento é executado pela função pdf(). Essa função recebe os mesmos argumentos da função png(). Note que assim como na função png(), também devemos informar a extensão ao nome do arquivo. Portanto, no fim do nome selcionado, devemos adicionar a extensão “.pdf”.
#Criando o arquivo .png que irá receber o gráfico
pdf(file="plot_modulo4_aula2_utilizando_terminal.pdf", width=10, height=10)
#Plotando objeto1 e objeto2 e editando alguns critérios dos eixos
plot(x=objeto1, y=objeto2, xlab="Sequência 1 (1 a 10)", ylab="Sequência 2 (10 a 100)", main="Gráfico simples", cex.lab=1.5, cex.axis=2, xlim=c(0,10), ylim=c(0,100), type="b", pch=17, cex=2, col="red", bty="n", axes=F)
#Criando o novo eixo-x contendo valores incrementados de 1 em 1
axis(side=1, at=seq(1,10,1), labels=seq(1,10,1), cex.axis=0.7)
#Criando o novo eixo-y contendo valores incrementados de 10 em 10
axis(side=2, at=seq(10,100,10), labels=seq(10,100,10), cex.axis=0.7)
#Resetando os prâmetros gráficos anteriores
dev.off()## png
## 2IMPORTANTE: Note que diferentemente da função png(), a unidade padrão para os valores informados aos argumento “width=” e “height=” é polegadas e não pixels. Por isso utilizamos os valores 10 e 10 para esses argumentos.
Verifique o arquivo “plot_modulo4_aula2_utilizando_terminal.pdf” presente na pasta do Dropbox para visualizar o output dos comandos anteriores.
À partir dos comandos acima, crie novos comandos alterando o nome dos arquivos salvos, os formatos e suas dimensões. Verifique o que é alterado por meio de cada modificação feita por você.
Criação de gráficos constantemente utilizados em rotinas de pesquisa, por meio do R
Nos exemplos anteriores, utilizamos dois objetos simples para aprentar os argumentos básicos para edição de plots no R. À partir de agora, iremos utilizar a base de dados “Modulo4_DataBase_CrohnD.txt” para a criação de alguns plots muito úteis em rotinas de pesquisa. Além disso, iremos abordar alguns tópicos de estatística descritiva não abordados na aula antrior devido a necessidade de criação de gráficos para potencializar a interpretação dos resultados.
Primeiramente, iremos inserir a base de dados “Modulo4_DataBase_CrohnD.txt” na variável crohn_data.
#Determinando diretório de trabalho
setwd("E:\\Doutorado\\curso_R\\Aula4")
#Carregando a base de dados para a variável
crohn_data<-read.table("Modulo4_DataBase_CrohnD.txt", h=T, sep="\t")
#Verificando a importação
head(crohn_data)## ID nrAdvE BMI height country sex age weight treat
## 1 19908 4 25.22 163 c1 F 47 67 placebo
## 2 19909 4 23.80 164 c1 F 53 64 d1
## 3 19910 1 23.05 164 c1 F 68 62 placebo
## 4 20908 1 25.71 165 c1 F 48 70 d2
## 5 20909 2 25.95 170 c1 F 67 75 placebo
## 6 20910 2 28.70 168 c1 F 54 81 d1Gráfico de pontos para representar correlação
A estimativa de correlação entre duas variáveis pode ser algo muito útil para se entender a relação entre as variáveis presentes em uma base de dados. No R, a correlação entre duas variáveis é calculada por meio da função cor.test(). É importante ressaltar que existem outras funções que também realizam testes de correlação. Contudo, a função cor.test() foi escolhida devido ao fato de executar também um teste de siginificância para a correlação. A função cor.test() recebe como argumentos “x=” e “y=”. Esses argumentos recebem os conjuntos de valores em que desejamos testar a correlação. Além disso, podemos selecionar as hipóteses alternativas para o teste por meio do argumento “alternative=” (“two.sided”, “greater”, ou “less”); o metodo para estimar a correlação por meio do argumento “method=” ( “pearson”, “kendall”, ou “spearman”), dentre outros argumentos. Contudo, no exemplo a seguir utilizaremos as definições padrão da função.
A seguir, vamos calcular a correlação entre o índice de massa corporal e altura, presentes nas coluna “BMI” e “heigth” da Data frame crohn_data. O resultado da correlação será inserido na variável “cor.result”.
#Calculando a correlação entre BMI e heigth
cor.result<-cor.test(x=crohn_data$BMI, crohn_data$height)
#Exibindo o resultado de cor.result
cor.result##
## Pearson's product-moment correlation
##
## data: crohn_data$BMI and crohn_data$height
## t = -1.1928, df = 115, p-value = 0.2354
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.28633347 0.07244009
## sample estimates:
## cor
## -0.1105466#Verificando a estrutura de cor.result
str(cor.result)## List of 9
## $ statistic : Named num -1.19
## ..- attr(*, "names")= chr "t"
## $ parameter : Named int 115
## ..- attr(*, "names")= chr "df"
## $ p.value : num 0.235
## $ estimate : Named num -0.111
## ..- attr(*, "names")= chr "cor"
## $ null.value : Named num 0
## ..- attr(*, "names")= chr "correlation"
## $ alternative: chr "two.sided"
## $ method : chr "Pearson's product-moment correlation"
## $ data.name : chr "crohn_data$BMI and crohn_data$height"
## $ conf.int : atomic [1:2] -0.2863 0.0724
## ..- attr(*, "conf.level")= num 0.95
## - attr(*, "class")= chr "htest"Note que o output da função cor.test() apresenta um valor de p para a significância da correlação, a hipotese alternativa para o teste de significância realizado, o intervalo de confiância e o valor da correlação. Observamos uma pequena correlação negativa (cor=-0.1105466) e não siginificativa (p-value=0.2354) entre as variáveis. Uma forma de representar graficamente essa correlação é por meio da utilização do gráfico de pontos entre as duas variáveis.
#Plotando BMI e heigth
plot(x=crohn_data$BMI, y=crohn_data$heigth, xlab="BMI", ylab="Heigth", cex.lab=1.5, cex.axis=2, pch=16, cex=2, bty="n")
Existem algumas edições que podem ser feitas no gráfico apresentado acima, para que o mesmo fique ainda mais informativo. Podemos adicionar uma “trend line”, indicando a tendência da correlação entre os dados. Além disso, podemos plotar um texto contendo os calores da correlação e o valor de p. Essas informações podem ser adicionadas pelas funções abline() e text().
A função abline() permite que a linha de tendência entre os dados seja plotada. Para isso, devemos informar como argumento principal o resultado de um modelo de regressão linear entre as variáveis. Nesse curso, não iremos entrar em detalhes quanto a interpretação do output de modelos lineares. Porém, para criarmos a linha de tendência, iremos realizar a regressão linear entre “BMI” e “height”, por meio da função lm(). O resultado dessa função será o primeiro argumento para a função abline(). Além disso, também iremos informar a cor e a espessura da linha que iremos plotar por meio dos argumentos “col=” e “lwd=”.
A função text(), recebe como argumentos básicos as posições no eixo-x e eixo-y onde será plotado o texto informado. Essas posições são informadas por meio dos argumentos “x=” e “y=”. Além disso, devemos informar o texto a ser plotado. O texto deve ser informado para o argumento “labels=”. Nesse exemplo, iremos extrair o valor da correlação diretamente da lista contendo os resultados do comando cor.test(), veja o output da função str() para mais detalhes. De modo a combinar textos com os valores extraídos do output armazendado na variável “cor.result”, iremos utilizar a função paste(). A função paste recebe como argumentos combinações de textos e objetos que serão combinados por um padrão que deve ser informado para o argumento “sep=”.
#Plotando BMI e heigth
plot(x=crohn_data$BMI, y=crohn_data$heigth, xlab="BMI", ylab="Heigth", cex.lab=1.5, cex.axis=2, pch=16, cex=2, bty="n")
#Plotando a linha de tendência
abline(lm(crohn_data$BMI~crohn_data$height), col="red", lwd=2)
#Plotando o texto com os valores de correlação e p-value
text(x=20, y=42, labels=paste("cor= ",round(cor.result$estimate,3), "; ", "p-value= ", round(cor.result$p.value,3), sep=""))
Note o gráfico gerado pelos comandos acima. Perceba os novos elementos. Nesse ponto do curso, atingimos um nível de combinação entre funções, estruturas de recuparação de dados em variáveis e programação que já é bastante complexo. Repare que a função round() foi utilizada para tornar os números mais apresentáveis por meio do arredondamento para 3 casas decimais.
DICA: Note que a linha de tendência reflete bem o valor de correlação. A inclinação da linha é bem leve e possuí um coeficiente de inclinação (o sentido do decaimento da curva) negativo. Isso é condizente com o valor de correlação obtido (-0.111).
Adicionalmente às informações que foram inseridas nos comandos anteriores, podemos alterar os caracteres utilizados no plot, em função de alguma das classes de nossa base de dados. Por exemplo, vamos alterar o plot caracter e a cor do gráfico em função do sexo dos participantes. Iremos fazer isso adicionando mais opções aos argumentos “pch=” e “col=”. Porém, deveremos utilizar uma estrutura de conexão entre as opções e as classes de nossa base de dados. Isso será feito com a utilização do “” logo após as opções informadas para cada argumento. Vamos determinar que homens serão representados por triângulos azuis e mulheres por círculos vermelhos. Além disso, um novo comando, que insere um novo parâmetro gráfico será apresentado. Nesse gráfico iremos inserir uma legenda, por meio da função legend(). Essa função deve ser utilizada logo após a criação do gráfico. A função legend() recebe como argumentos “x=” e “y=”, que irão informar as coordenadas para que a legenda seja plotada. Contudo, “x=” e “y=” podem ser substituídos pelas palavras: “topright”, “topleft”, “bottomright” e “bottomleft”. Essas palavras irão determinar que a legenda seja plotada no canto superior esquerdo, canto superior direito, canto inferior esquerdo e canto inferior direito, respectivamente. Além disso, devemos informar as cores e qual caracter desejamos utilizar para representar as cores na legenda. Essas escolhas são feitas por meio dos argumentos já conhecidos: “col=” e “pch=”, respectivamente. Por fim, devemos informar o texto que será adicionado para cada caracter da legenda. Isso é feito por meio do uso do argumento “legend=”. Adicionalmente, existem outros argumentos que podem ser utilizados, como para eliminar as bordas em torno da legenda, adicionar/alterar título, mudar a espessura da bordar em torno da legenda, preenchimento, etc.
#Plotando BMI e heigth
plot(x=crohn_data$BMI, y=crohn_data$heigth, xlab="BMI", ylab="Heigth", cex.lab=1.5, cex.axis=2, pch=c(16,17)[crohn_data$sex], col=c("red","blue")[crohn_data$sex], cex=2, bty="n")
#Plotando a linha de tendência
abline(lm(crohn_data$BMI~crohn_data$height), col="red", lwd=2)
#Plotando o texto com os valores de correlação e p-value
text(x=20, y=42, labels=paste("cor= ",round(cor.result$estimate,3), "; ", "p-value= ", round(cor.result$p.value,3), sep=""))
#Plotando a legenda no gráfico no canto superos direito
legend("topright", pch=c(16,17), col=c("red","blue"), legend=c("Mulheres", "Homens"))
Note como “[crohn_data$sex]” foi utilizado logo após os valores informados aos argumentos “pch=” e “col=”. Essa estrutura permite que seja feita a conexão entre os valroes informados e as classes. Contudo, sempre é necessário voltar ao output da função str(), uma vez que o primeiro valor informado ao argumento, será conectado à primeira classe presente no vetor contendo os grupos. Por exemplo, se avaliarmos o output de str(crohn_data), veremos que a primeira classe da coluna “sex” é “F”, de female. Portanto, a cor “red”, informada para o argumento “col=” será associada às mulheres, assim consequentemente para as outras classes, no mesmo padrão.
DETALHE: Os gráficos de pontos são constatemente chamados de “scatter plots”. Geralmente esses gráficos são utilizados para representar o padrão de relação entre quaisquer duas variáveis.
Testes de normalidade e histogramas
A normalidade de um dado reflete a distribuição do mesmo. Mais especificamente, se esse dados está dentro de uma distribuição normal. Como o objetivo desse curso não é discutir conteúdo estatístico, mais informações sobre a distribuição normal podem ser encontradas aqui: http://www.portalaction.com.br/probabilidades/62-distribuicao-normal
Contudo, é importante saber que identificar se um dado possui ou não uma distribuição normal é um pressuposto para vários testes estatísticos, como por exemplo, o teste de Bartlett e a ANOVA.
Existe uma série de testes para avaliar se um número está dentro de uma distribuição normal ou não. Aqui, iremos utilizar a função shapiro.test(). Essa função executa o teste de normalidade de Shapiro-Wilk para avaliar se o dado está em uma distribuição normal.Além dos testes, é sempre recomendado que seja feita uma inspeção visual da distribuição do dado. Geralmente, essa insperção é feita por meio de histogramas. Os histogramas são criado no R, por meio da função hist(). Abaixo, iremos comparar dois tipos de dados. Primeiramente iremos simular uma distribuição normal com média 5 e desvio padrão 3, contendo 100 valores por meio da função rnorm(). Além disso, também iremos avaliar a normalidade dos valores contidos em “BMI” na base de dados crohn_data. Em seguida iremos plotar os histogramas de ambos os valores, lado a lado.
#Criando um grupo de valores com distribuição normal, com média igual a 5 e désvio-padrão igua a 3
norm.data<-rnorm(100, mean = 5, sd = 3)
#Avaliando a normalidade do dado norm.data por meio da função shapiro.test
shapiro.test(rnorm(100, mean = 5, sd = 3))##
## Shapiro-Wilk normality test
##
## data: rnorm(100, mean = 5, sd = 3)
## W = 0.97506, p-value = 0.05459#Avaliando a normalidade do dado BMI por meio da função shapiro.test
shapiro.test(crohn_data$BMI)##
## Shapiro-Wilk normality test
##
## data: crohn_data$BMI
## W = 0.92364, p-value = 5.061e-06#Plotando ambos os histogramas, lado a lado
#Determinando os parâmetros gráficos
par(mfrow=c(1,2))
#Plotando o histograma para norm.data
hist(x=norm.data, ylab="Frequência", xlab="Observações", main="Histograma de norm.data", cex.axis=2, cex.lab=1.5)
#Plotando o histograma para BMI
hist(x=crohn_data$BMI, ylab="Frequência", xlab="Observações", main="Histograma de BMI", cex.axis=2, cex.lab=1.5)
#Resetando parâmetros gráficos
dev.off()## null device
## 1Para interpretarmos os resultados do teste de Shapiro-Wilk, precisamos saber que a H0 desse teste considera que o dado possui uma distribuição normal. Portanto, no caso da presença de valores de p que estejam abaixo do valor de alfa pré-determinado (geralmente, p<0.05), rejeitamos H0 e consideramos que os dados não seguem uma distribuição normal.
Note que para o objeto norm.data, o valor de p para o teste de Shapiro-Wilk foi 0.8777, assim, não rejeitando H0 e considerando o dado como possuindo uma distribuição normal. Isso já era esperado, uma vez que criamos o dado com essa característica. Já para os valores em BMI, obtivemos um valor de p igual a 5.061e-06, o que é muito menor que 0.05. Dessa forma, rejeitamos H0. Assim, considerando que o dado não segue uma distribuição normal.
Adicionalmente, quando comparamos os histogramas de ambos os dados podemos notar que o histograma de norm.data possui um padrão que se assemelha muito mais ao formato de sino, característico da distribuição normal. Já para o histograma de BMI, esse padrão de sino não é observado. Além disso, somos capazes de perceber uma assimetria entre as extremidades do gráfico, algo que é bastante decisivo para determinarmos que o dado não possui distribuição normal. O que inclusive, pode ser uma das causas para um valor de p tão baixo obtido no teste de Shapiro-Wilk.
DICA: Geralmente é aconselhado a aplicação de mais de um teste de normalidade para avaliar a distribuição do dado. Outros exemplos de teste de normalidade presentes no R são o teste de Anderson-Darling e o teste de Kolmogorov-Smirnov. Esses testes podem ser executados pelas funções ad.test() e ks.test(), respectivamente. A utilização da função ad.test() depende da instalação do pacote “nortest”.
Gráficos de barra utilizando o R
Um tipo de gráfico amplamente utilizado para a representação de valores entre diferentes classes presentes em uma base de dados, é o gráfico de barras. No R, o gráfico de barras pode ser gerado por meio da função barplot(). O funcionamento da função barplot() é um pouco distinto dos gráficos apresentados anteriormente. Ao invés de apresentarmos os argumentos que serão plotados nos eixos x e y, devemos informar à função a altura de cada uma das barras que serão plotadas no gráfico. Esses valores devem ser informados ao argumento “heigth=”. Além disso, um outro argumento pode ser informado à função barplot() para determinar o nome presente abaixo de cada barra. Esse argumento é o “names.arg=”. Outros argumentos podem ser informados à função barplot(), a maioria deles trata-se de argumentos já apresentados nessa aula e que realizam alterações na estrutura e estilo do plot. Porém, existem alguns argumentos específicos da função barplot() que podem ser alterado. Para mais detalhes, utilize o comando “?barplot()”.
Como exemplo da utilização da função barplot(), vamos plotar um gráfico de barras contendo o número de participantes presentes em cada classe da coluna “treat”, da data base crohn_data. Para a obtenção da quantidade de indivíduos presentes em cada tratamento, iremos utilizar a combinação das funções table() e as.data.frame(), já apresentadas na aula anterior.
#Criando uma Data frame com o número de indivíduos para cada tratamento
data.treat<-as.data.frame(table(crohn_data$treat))
#Visualizando a Data frame criada
data.treat## Var1 Freq
## 1 d1 39
## 2 d2 39
## 3 placebo 39#Plotando o gráfico de barras para os valores obtidos
barplot(height = data.treat$Freq, names.arg = c("Droga 1", "Droga 2", "Placebo"), ylab="Número de indivíduos", xlab="Tipo de tratamento", cex.axis=2, cex.lab=1.5,bty="n", ylim=c(0,50))
Note que assim como já observado anteriormente, todos os grupos possuem o mesmo número de indivíduos (39). Além disso, note que na função barplot(), optei por utilizar o argumento “ylim=”. Para esse argumento utilizei como valor mínimo 0 e valor máximo 50. Isso foi devido ao fato do valor máximo presente no argumento “heigth=” ser 39 e a incrementação do eixo-y, pré-determianda pelo R, ser de 10 em 10. Como na forma padrão o R utiliza o maior valor como o limite superior do eixo, não é possível incrementar até 40, já que o valor máximo é 39. Portanto, na ausência do argumento ylim=c(0,50), o eixo-y teria sido incrementado somente até 30. Tente executar o comando acima retirando o argumento ylim=c(0,50) para visualizar o resultado. Além disso, tente mudar a cor das barras do plot.
Uma forma intressante de se criar gráficos de barras é plotando as barras para observações de mais de um grupo, lado a lado. Por exemplo, vamos plotar o número de particpantes de cada tratamento em função do sexo. Para isso, vamos usar o output da função table() para as duas variáveis, sexo e tratamento. Isso será informado para o argumento “heigth=”. Além disso, vamos colorir as barras utilizando o mesmo padrão de cores utilziados no gráfico de pontos.
#Obtendo a tabela com os valores por cada combinação entre tratamento e sexo
counts <- table(crohn_data$sex, crohn_data$treat)
#Visualizando o conteúdo de counts
counts##
## d1 d2 placebo
## F 34 34 32
## M 5 5 7#Plotando o gráfico de barras para os valores obtidos
barplot(height= counts, names.arg = c("Droga 1", "Droga 2", "Placebo"), col=c("red", "blue"), ylab="Número de indivíduos", xlab="Tipo de tratamento", cex.axis=2, cex.lab=1.5,bty="n", ylim=c(0,50), beside=T)
#Plotando a legenda no gráfico no canto superos direito
legend("topright", pch=c(15,15), col=c("red","blue"), legend=c("Mulheres", "Homens"))
Note que para gerar esse tipo de gráfico, fizemos duas alterações. Primeiro, foi informado a tabela contendo a quantidade de indivíduos para cada combinação entre sexo e tratamento, obtido por meio do comando table(). Além disso, um novo argumento foi inserido em barplot. Esse argumento foi o “beside=”. Quando inserimos uma resposta lógica igual a verdadeiro (“TRUE”) para esse argumento, plotamos as diferentes colunas presentes em counts lado a lado no plot.
Gráficos de pizza utilizando o R
O gráfico de pizza é constantemente utilizado para representar as proporções entre diferentes classes em uma base de dados. No R, os gráficos de pizza são criados por meio da função pie(). Essa função recebe como argumentos principais “x=” e “labels=”. Em “x=”, devemos informar as proporções de cada fatia do gráfico de pizza. Para “labels=”, são informados os nomes de cada fatia. Além desses argumentos, podemos informar o argumento “col=”. Com esse argumento podemos determinar as cores que serão utilizadas em cada fatia. Assim como para as demais funções de plot apresentadas até aqui, outros argumentos específicos são possíveis de serem utilizados. Mais informações podem ser sempre encontradas por meio da leitura do manual da função. No caso da função pie(), o manual pode ser acessado utilizando o comando “?pie()”.
Como exemplo da criação de gráficos de pizza no R, iremos criar um gráfico de barras para as proporções entre os sexos na base de dados crohn_data. Para calcular as proporções, iremos utilizar a função table(), já apresentada na aula anterior. Além disso, iremos deixar os nomes de cada fatia do gráfico um pouco mais informativos por meio do uso da função paste(). Para a porção do plot referente aos homens, iremos utilizar azul. Para a porção do plot referente às mulheres, iremos utilizar vermelho.
#Criando uma Data frame com as proporções de cada sexo na base de dados
data.prop<-as.data.frame(table(crohn_data$sex))
#Verificando o output inserido em data.prop
data.prop## Var1 Freq
## 1 F 100
## 2 M 17#Criando as proporções para serem utilizadas em labels
labels.prop<-c(round(((data.prop[,2])/sum(data.prop[,2]))*100,2))
#Criando as porcentagens para utilizarmos no argumento labels
labels.comp<-paste(labels.prop, "%", sep="")
#Criando o gráfico de pizza para as proporções contidas em data.prop
pie(x=c(data.prop[1,2],data.prop[2,2]), labels=labels.comp, col=c("red", "blue"), main="Proporção entre homens e mulheres em crohn_data")
#Plotando a legenda no gráfico no canto superos direito
legend("topright", pch=15, col=c("red","blue"), legend=c("Mulheres", "Homens"))
Box-plots utilizando o R
Os box-plots são uma classe de gráficos muito úteis para a apresentação de dados em estatística descritiva. Com o box-plot somos capazes de representar uma série e elementos que descrevem bem o comportamento de nossos dados. No box-plot conseguimos retirar informações referentes à mediana, primeiro quatil, terceiro quartil, o limite superior, o limite inferior e presença de outliers. Mais informações sobre a interpretação de box-plots podem ser encontradas aqui: https://pt.wikipedia.org/wiki/Diagrama_de_caixa
No R, os box-plots são criados por meio da função boxplot(). Geralmente, o box-plot é utilizado para a comparação de diversos grupos em um mesmo gráficos. Aqui, iremos comparar o número de eventos adversos para cada grupo de tratamento. Para isso, a função boxplot() recebe como argumento uma fórmula que representa a relação entre o conjunto de dados de número de eventos adversos e o conjunto de dados que determina cada observação para um dos grupos. Essa relação é feita por meio da tulização do caracter já descrito previamente, “~”. Contudo, caso alguém deseje criar um box-plot contendo somente a representação de um grupo ou de um conjunto de dados sem grupos, isso pode ser feito de uma meneira muito simples. Para isso, devemos informar somente o conjunto de dados para o argumento “x=”, na função boxplot().
#Criando o box-plot para o número de reações adversas para cada grupo de tratamento
boxplot(crohn_data$nrAdvE ~ crohn_data$treat, names=c("Droga 1", "Droga 2", "Placebo"),ylab="Número de reações adversas", xlab="Observações", cex.axis=2, cex.lab=1.5)
Observe o output criado. Podmeos notar nos box-plots que a mediana (linha preta dentro da caixa) do número de rações adversas para a droga 1 é menor do que a mediana para a droga 1 e o placedo. No entanto, entre a droga 2 e o placebo, as medianas são iguais. A altura das caixas representa os valores do primeiro quartil (representado pela parte inferior da caixa) e o terceiro quartil (representado pela parte inferior da caixa). As linhas e barras acima das caixas representam os respectivos limites superiores para um intervalo de confiança de 95%. Repare que não há linhas abaixo, pois o primeiro quartil coincide com o limite inferior do intervalo de confiânça. Além disso, conseguimos notar alguns pontos no box-plot. Esses pontos são os chamados outliers, valores que estão fora do intervalo de confiânça de 95%. Da mesma forma que fizemos para o plot anterior, podemos colorir os box-plots em função das classes. Tente executar essa ação.
Considerações finais
Nesse módulo foi apresentado uma série de funções relacionadas à criação de gráficos no ambiente R. Existem uma infinidade de outros gráficos e edições possíveis de serem feitos no R. Portanto, os próximos passos para o aprendizado de novos comandos devem ser guiados pela demanda de cada um de vocês. Aqui, cobrimos uma boa parte dos comandos básicos, o que fará o processo de aprendizado desses novos comandos algo bem mais fácil. Além disso, existem outros pacotes criados somente para a produção de gráficos mais complexos ou com uma estetíca mais aprimorada. Alguns exemplos são o “ggplot2” e o “plotly”. Esses pacotes possuem os seus próprio grupos de funções e argumentos, assim como as específicidades de seus outputs. Cada um desses pacotes poderiam ser o assunto de cursos específicos, dado a grande variedade de funções.
Um site muito útil para aprender sobre diferentes plots e sobre o uso de outras ferramentas para criação de gráficos no R é o R-graph gallery. Nesse site, podemos encontrar diversos tipos de plots explicados e com os códigos utilizados para criá-los. Procure acesar o site, copiar e colar os comandos em seu terminal e ir isolando cada um dos comandos contidos na comando maior. Com isso, você conseguirá notar o que cada uma das informações passadas como argumentos alteram no gráfico final. Como dito no início do curso, o processo de aprendizagem envolvendo programação demanda treino e memorização. Para isso, devemos sempre praticar o uso das novas ferramentas aprendidas.