Introdução ao R e aplicações em Genética
Esse material foi elaborado para o Introdução ao R e aplicações em Genética.
Sugerimos que, antes de iniciar a prática aqui descrita, siga este tutorial para instalação do R e do RStudio.
Familiarização com a interface do RStudio
Abrindo o RStudio você verá:
A interface é separada em quatro janelas com principais funções:
- Edição de código
- Ambiente de trabalho e histórico
- Console
- Arquivos, gráficos, pacotes e ajuda
Explore cada uma das janelas. São inúmeras funcionalidades para cada uma delas, veremos algumas delas no decorrer do curso.
Um primeiro script
A janela de edição de código (provavelmente localizada no canto superior esquerdo) será utilizada para escrever o seu código. Abra um novo script clicando no + no canto superior esquerdo e selecionando R script.
Vamos então iniciar os trabalhos com o tradicional Hello World. Digite no seu script:
cat("Hello world")## Hello worldAgora, selecione a linha e aperte o botão Run ou utilize Ctrl + enter.
Ao fazer isso o seu código será processado na janela Console, onde aparecerá em azul (se você estiver com as cores padrão do R) o código escrito e, logo em seguida, o resultado desejado. A linha somente não será processada no console se houver o símbolo # na frente. Agora, experimente colocar # na frente do código escrito. Novamente, selecione a linha e aperte Run.
# cat("Hello world")O símbolo # é utilizado para comentários no código. Esta é uma ótima prática de organização e ajuda a lembrar, posteriormente, o que você estava pensando quando escreveu o código. Também é essencial para que outras pessoas possam entendê-lo. Como no exemplo:
# Iniciando os trabalhos no R
cat("Hello world")## Hello worldImportante: sempre que quiser realizar alguma alteração, edite o seu script e não diretamente no console, pois tudo o que neste é escrito, não terá como ser salvo!
Para salvar seu script, você pode utilizar a aba Files localizada (como padrão) no canto direito inferior. Você pode procurar uma localização de sua preferência, criar uma nova pasta com o nome CursoR.
Dica:
- Evite colocar espaços e pontuações no nome das pastas e arquivos, isso pode dificultar o acesso via linha de comando no R. Por exemplo, ao invés de
Curso R, optamos porCursoR.
Depois, basta clicar no disquete localizado no cabeçalho do RStudio ou com Ctrl + s e selecionar o diretório CursoR criado. Scripts em R são salvos com a extensão .R.
Estabelecendo diretório de trabalho
Outra boa prática no R é deixar o script no mesmo diretório onde estão seus dados brutos (arquivos de entrada no script) e os dados processados (gráficos, tabelas, etc). Para isso, vamos fazer com que o R identifique o mesmo diretório em que você salvou o script como diretório de trabalho. Assim, ele entenderá que é dali que os dados serão obtidos e é para lá que também irão os resultados.
Você pode fazer isso utilizando as facilidades do RStudio, basta localizar o diretório CursoR pela aba Files, clicar em More e depois “Set as Working Directory”. Repare que irá aparecer no console algo como:
setwd("~/Documents/CursoR")Ou seja, você pode utilizar este mesmo comando para realizar esta ação. O resultado será nossa pasta de trabalho. Quando estiver perdido/a ou para ter certeza que o diretório de trabalho foi alterado utilize:
getwd()Facilitando a vida com Tab
Agora, imagine que você tem um diretório como ~/Documentos/mestrado/semestre1/disciplina_tal/aula_tal/dados_28174/analise_276182/resultados_161/. Não é fácil lembrar todo este caminho para escrever num comando setwd().
Além da facilidade da janela do RStudio, você também pode utilizar a tecla Tab para completar o caminho para você. Experimente buscando alguma pasta no seu computador. Basta começar a digitar o caminho e apertar Tab, ele irá completar o nome para você! Se você tiver mais do que um arquivo com aquele início de nome, aperte duas vezes o Tab, ele mostrará todas as opções.
O Tab funciona não só para indicar caminhos, mas também para comandos e nomes de objetos. É muito comum cometermos erros de digitação no código. Utilizar o Tab reduzirá significativamente esses erros.
Operações básicas
Vamos então à linguagem!
O R pode funcionar como uma simples calculadora, que utiliza a mesma sintaxe que outros programas (como o excel):
1+1.3 #Decimal definido com "."
2*3
2^3
4/2
sqrt(4) #raíz quadrada
log(100, base = 10) #Logaritmo na base 10
log(100) #Logaritmo com base neperianaAgora, utilize as operações básicas para solucionar expressão abaixo. Lembre-se de utilizar parênteses () para estabelecer prioridades nas operações.
\((\frac{13+2+1.5}{3})+ log_{4}96\)
Resultado esperado:
## [1] 8.792481Repare que, se posicionar o parênteses de forma incorreta, o código não resultará em nenhuma mensagem de erro, pois este é um erro que chamamos de erro lógico, ou seja, o código roda, mas não faz o que você gostaria que ele fizesse. Esse é o tipo de erro mais difícil de ser consertado. Os erros que produzem uma mensagem, seja um aviso (warning) ou um erro (error) são chamados de erros de sintaxe. Nesses casos, o R retornará uma mensagem para te ajudar a corrigí-los. Os warnings não comprometem o funcionamento do código, mas chamam a atenção para algum ponto; já os errors precisam necessariamente ser corrigidos para que o código rode.
Exemplo de error:
((13+2+1,5)/3) + log(96, base = 4)Você pode também esquecer de fechar algum parênteses, ou aspas, ou colchetes, ou chaves, nesses casos, o R ficará aguardando o comando para fechar o bloco de código sinalizando com um +:
((13+2+1.5)/3 + log(96, base = 4)Se acontecer, vá até o console e aperte ESC, que o bloco será finalizado para que você possa corrigí-lo.
Os comandos log e sqrt são duas de muitas outras funções básicas que o R possui. Funções são conjuntos de instruções organizadas para realizar uma tarefa. Para todas elas, o R possui uma descrição para auxiliar no seu uso. Para acessar essa ajuda use:
?logE será aberta a descrição da função na janela Help do RStudio.
Se a descrição do próprio R não for suficiente para você entender como funciona a função, busque no google (de preferência em inglês). Existem diversos sites e fóruns com informações didáticas das funções do R.
Operações com vetores
Os vetores são as estruturas mais simples trabalhadas no R. Construímos um vetor com uma sequencia numérica usando:
c(1,3,2,5,2)## [1] 1 3 2 5 2MUITA ATENÇÃO: O c é a função do R (Combine Values into a Vector or List) com a qual construímos um vetor!
Utilizamos o simbolo : para criar sequências de números inteiros, como:
1:10## [1] 1 2 3 4 5 6 7 8 9 10Podemos utilizar outras funções para gerar sequências, como:
seq(from=0, to=100, by=5)## [1] 0 5 10 15 20 25 30 35 40 45 50 55 60 65 70 75 80 85 90
## [20] 95 100# ou
seq(0,100,5) # Se você já souber a ordem dos argumentos da função## [1] 0 5 10 15 20 25 30 35 40 45 50 55 60 65 70 75 80 85 90
## [20] 95 100- Crie uma sequencia utilizando a função
seqque varie de 4 a 30, com intervalos de 3 em 3.
## [1] 4 7 10 13 16 19 22 25 28A função rep gera sequências com números repetidos:
rep(3:5, 2)## [1] 3 4 5 3 4 5Podemos realizar operações utilizando esses vetores:
c(1,4,3,2)*2
c(4,2,1,5)+c(5,2,6,1)
c(4,2,1,5)*c(5,2,6,1)Repare que já esta ficando cansativo digitar os mesmos números repetidamente, vamos resolver isso criando objetos para armazenar nossos vetores e muito mais.
Criando objetos
O armazenamento de informações em objetos e a possível manipulação desses faz do R uma linguagem orientada por objetos. Para criar um objeto basta atribuir valores para as variáveis, como a seguir:
x = c(30.1,30.4,40,30.2,30.6,40.1)
# ou
x <- c(30.1,30.4,40,30.2,30.6,40.1)
y = c(0.26,0.3,0.36,0.24,0.27,0.35)Os mais antigos costumam usar o sinal <-, mas tem a mesma função de =. Há quem prefira usar o <- como atribuição em objetos e = apenas para definir os argumentos dentro de funções. Organize-se da forma como preferir.
Para acessar os valores dentro do objeto basta:
x## [1] 30.1 30.4 40.0 30.2 30.6 40.1A linguagem é sensível à letras maiúsculas e minúsculas. Logo, x é diferente de X:
XO objeto X não foi criado.
O nome dos objetos é uma escolha pessoal, a sugestão é tentar manter um padrão para melhor organização. Alguns nomes não podem ser usados por estabelecerem papéis fixos no R, são eles:
- NA - Not available, simboliza dados faltantes
- NaN - Not a number, simboliza indefinições matemáticas
- Inf - Infinite, conceito matemático
- NULL - Null object, simboliza ausência de informação
Podemos então realizar as operações com o objeto criado:
x+2## [1] 32.1 32.4 42.0 32.2 32.6 42.1x*2## [1] 60.2 60.8 80.0 60.4 61.2 80.2Para realizar a operação o R alinha os dois vetores e realiza a operação elemento à elemento. Observe:
x + y## [1] 30.36 30.70 40.36 30.44 30.87 40.45x*y## [1] 7.826 9.120 14.400 7.248 8.262 14.035Se os vetores tiverem tamanhos diferentes, ele irá repetir o menor para realizar a operação elemento a elemento com todos do maior.
x*2
x*c(1,2)Se caso o menor vetor não for múltiplo do maior, obteremos um aviso:
x*c(1,2,3,4)## Warning in x * c(1, 2, 3, 4): longer object length is not a multiple of shorter
## object length## [1] 30.1 60.8 120.0 120.8 30.6 80.2Repare que o warning não compromente o funcionamento do código, ele só dá uma dica de que algo pode não estar da forma como você gostaria.
Podemos também armazenar a operação em outro objeto:
z <- (x+y)/2
zPodemos também aplicar algumas funções, como exemplo:
sum(z) # soma dos valores de z## [1] 101.59mean(z) # média ## [1] 16.93167var(z) # variância## [1] 6.427507Indexação
Acessamos somente o 3º valor do vetor criado com []:
z[3]Também podemos acessar o número da posição 2 a 4 com:
z[2:4]## [1] 15.35 20.18 15.22Para obter informações do vetor criado utilize:
str(z)## num [1:6] 15.2 15.3 20.2 15.2 15.4 ...A função str nos diz sobre a estrutura do vetor, que se trata de um vetor numérico com 6 elementos.
Os vetores também podem receber outras categorias como caracteres:
clone <- c("GRA02", "URO01", "URO03", "GRA02", "GRA01", "URO01")Outra classe são os fatores, esses podem ser um pouco complexos de lidar.
De forma geral, fatores são valores categorizados por levels, como exemplo, se transformarmos nosso vetor de caracteres clone em fator, serão atribuidos níveis para cada uma das palavras:
clone_fator <- as.factor(clone)
str(clone_fator)## Factor w/ 4 levels "GRA01","GRA02",..: 2 3 4 2 1 3levels(clone_fator)## [1] "GRA01" "GRA02" "URO01" "URO03"Dessa forma, teremos apenas 4 níveis para um vetor com 6 elementos, já que as palavras “GRA02” e “URO01” se repetem. Podemos obter o número de elementos do vetor ou o seu comprimento com:
length(clone_fator)## [1] 6Também há vetores lógicos, que recebem valores de verdadeiro ou falso:
logico <- x > 40
logico # Os elementos são maiores que 40?## [1] FALSE FALSE FALSE FALSE FALSE TRUECom ele podemos, por exemplo, identificar quais são as posições dos elementos maiores que 40:
which(logico) # Obtendo as posiçoes dos elementos TRUE## [1] 6x[which(logico)] # Obtendo os números maiores que 40 do vetor x pela posição## [1] 40.1# ou
x[which(x > 40)]## [1] 40.1Também podemos localizar elementos específicos com:
clone %in% c("URO03", "GRA02")## [1] TRUE FALSE TRUE TRUE FALSE FALSETambém podem ser úteis as funções any e all. Procure sobre elas.
Encontre mais sobre outros operadores lógicos, como o > utilizado, neste link.
Warning1
Faça uma sequência numérica, contendo 10 valores inteiros, e salve em um objeto chamado “a”.
(a <- 1:10)## [1] 1 2 3 4 5 6 7 8 9 10Crie outra sequência, utilizando números decimais e qualquer operação matemática, de tal forma que seus valores sejam idênticos ao objeto “a”.
b <- seq(from = 0.1, to = 1, 0.1)
(b <- b*10)## [1] 1 2 3 4 5 6 7 8 9 10Os dois vetores parecem iguais, não?
Então, utilizando um operador lógico, vamos verificar o objeto “b” é igual ao objeto “a”.
a==b## [1] TRUE TRUE FALSE TRUE TRUE TRUE FALSE TRUE TRUE TRUEAlguns valores não são iguais. Como isso é possivel?
a==round(b)## [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUEWarning2
Não é possível misturar diferentes classes dentro de um mesmo vetor. Ao tentar fazer isso, repare que o R irá tentar igualar para uma única classe:
errado <- c(TRUE, "vish", 1)
errado## [1] "TRUE" "vish" "1"No caso, todos os elementos foram transformados em caracteres.
Algumas Dicas:
- Cuidado com a prioridade das operações, na dúvida, sempre acrescente parenteses conforme seu interesse de prioridade.
- Lembre-se que, se esquecer de fechar algum
(ou[ou", o console do R ficará esperando você fechar indicando um+. Nada será processado até que você digite diretamente no console um)ou aperte ESC. - Cuidado para não sobrepor objetos já criados criando outros com o mesmo nome. Use, por exemplo: altura1, altura2.
- Mantenha no seu script .R somente os comandos que funcionaram e, de preferência, adicione comentários. Você pode, por exemplo, comentar dificuldades encontradas, para que você não cometa os mesmos erros mais tarde.
Se estiver adiantada/o em relação aos colegas, você já pode fazer os exercícios da Sessão 1, se não, faça-os em outro momento e nos envie dúvidas.
Matrizes
As matrizes são outra classe de objetos muito utilizadas no R, com elas podemos realizar operações de maior escala de forma automatizada.
Por serem usadas em operações, normalmente armazenamos nelas elementos numéricos. Para criar uma matriz, determinamos uma sequência de números e indicamos o número de linhas e colunas da matriz:
X <- matrix(1:12, nrow = 6, ncol = 2)
X## [,1] [,2]
## [1,] 1 7
## [2,] 2 8
## [3,] 3 9
## [4,] 4 10
## [5,] 5 11
## [6,] 6 12Podemos também utilizar sequências já armazenadas em vetores para gerar uma matriz, desde que eles sejam numéricos:
W <- matrix(c(x,y), nrow = 6, ncol =2)
W## [,1] [,2]
## [1,] 30.1 0.26
## [2,] 30.4 0.30
## [3,] 40.0 0.36
## [4,] 30.2 0.24
## [5,] 30.6 0.27
## [6,] 40.1 0.35Com elas podemos realizar operações matriciais:
X*2## [,1] [,2]
## [1,] 2 14
## [2,] 4 16
## [3,] 6 18
## [4,] 8 20
## [5,] 10 22
## [6,] 12 24X*X ## [,1] [,2]
## [1,] 1 49
## [2,] 4 64
## [3,] 9 81
## [4,] 16 100
## [5,] 25 121
## [6,] 36 144X%*%t(X) # Multiplicação matricial## [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 50 58 66 74 82 90
## [2,] 58 68 78 88 98 108
## [3,] 66 78 90 102 114 126
## [4,] 74 88 102 116 130 144
## [5,] 82 98 114 130 146 162
## [6,] 90 108 126 144 162 180Utilizar essas operações exige conhecimento de álgebra de matrizes, se quiser se aprofundar a respeito, o livro Linear Models in Statistics, Rencher (2008) possui um boa revisão à respeito. Você também pode explorar a sintaxe do R para essas operações neste link.
Acessamos os números internos à matriz dando as coordenadas [linha,coluna], como no exemplo:
W[4,2] # Número posicionado na linha 4 e coluna 2## [1] 0.24As vezes pode ser informativo dar nomes às colunas e às linhas da matriz, fazemos isso com:
colnames(W) <- c("altura", "diametro")
rownames(W) <- clone
W## altura diametro
## GRA02 30.1 0.26
## URO01 30.4 0.30
## URO03 40.0 0.36
## GRA02 30.2 0.24
## GRA01 30.6 0.27
## URO01 40.1 0.35Essas funções colnames e rownames também funcionam nos data.frames.
Data.frames
Diferente das matrizes, não realizamos operações com os data.frames, mas eles permitem a união de vetores com classes diferentes. Os data frames são semelhantes a tabelas geradas em outros programas, como o excel.
Os data frames são combinações de vetores de mesmo comprimento. Todos os que criamos até agora tem tamanho 6, verifique.
Podemos assim combiná-los em colunas de um único data.frame:
campo1 <- data.frame("clone" = clone, # Antes do sinal de "="
"altura" = x, # estabelecemos os nomes
"diametro" = y, # das colunas
"idade" = rep(3:5, 2),
"corte"= logico)
campo1## clone altura diametro idade corte
## 1 GRA02 30.1 0.26 3 FALSE
## 2 URO01 30.4 0.30 4 FALSE
## 3 URO03 40.0 0.36 5 FALSE
## 4 GRA02 30.2 0.24 3 FALSE
## 5 GRA01 30.6 0.27 4 FALSE
## 6 URO01 40.1 0.35 5 TRUEPodemos acessar cada uma das colunas com:
campo1$idade## [1] 3 4 5 3 4 5Ou também com:
campo1[,4] ## [1] 3 4 5 3 4 5Aqui, o número dentro dos colchetes se refere à coluna, por ser o segundo elemento (separado por vírgula). O primeiro elemento se refere à linha. Como deixamos o primeiro elemento vazio, estaremos nos referindo a todas as linhas para aquela coluna.
Dessa forma, se quisermos obter um conteúdo específico podemos dar as coordenadas com [linha,coluna]:
campo1[1,2] ## [1] 30.1- Obtenha o diâmetro do clone “URO03”.
## [1] 0.36Mesmo se tratando de um data frame, podemos realizar operações com os vetores numéricos que a compõe.
- Com o diâmetro e a altura das árvores, calcule o volume conforme a fórmula a seguir e armazene em um objeto
volume:
\(3.14*(diametro/2)^2*altura\)
## [1] 1.597287 2.147760 4.069440 1.365523 1.751131 3.856116Agora, vamos adicionar o vetor calculado com o volume ao nosso data frame. Para isso use a função cbind.
campo1 <- cbind(campo1, volume)
campo1## clone altura diametro idade corte volume
## 1 GRA02 30.1 0.26 3 FALSE 1.597287
## 2 URO01 30.4 0.30 4 FALSE 2.147760
## 3 URO03 40.0 0.36 5 FALSE 4.069440
## 4 GRA02 30.2 0.24 3 FALSE 1.365523
## 5 GRA01 30.6 0.27 4 FALSE 1.751131
## 6 URO01 40.1 0.35 5 TRUE 3.856116str(campo1)## 'data.frame': 6 obs. of 6 variables:
## $ clone : chr "GRA02" "URO01" "URO03" "GRA02" ...
## $ altura : num 30.1 30.4 40 30.2 30.6 40.1
## $ diametro: num 0.26 0.3 0.36 0.24 0.27 0.35
## $ idade : int 3 4 5 3 4 5
## $ corte : logi FALSE FALSE FALSE FALSE FALSE TRUE
## $ volume : num 1.6 2.15 4.07 1.37 1.75 ...Algumas dicas:
Lembre-se que, para construir matrizes e data frames, as colunas devem apresentar o mesmo número de elementos.
Caso não saiba o operador ou a função que deve ser utilizada, busque no google. Por exemplo, se está em dúvida em como calcular o desvio padrão, busque por “desvio padrão R” ou, melhor ainda, “standard deviation R”. Logo nas primeiras páginas você obterá respostas. A comunidade do R é bastante ativa e grande parte das suas perguntas sobre ele já foram respondidas em algum lugar da web.
Não esqueça que tudo o que fizer no R precisa ser explicitamente indicado, como uma multiplicação 4ac com
4*a*c. Para gerar um vetor 1,3,2,6 é necessário:c(1,3,2,6).
Paramos aqui no primeiro dia!
Obtenha o arquivo RData com todos os objetos criados até então no tutorial aqui.
Obtenha o arquivo Rscript com os comandos rodados até então no tutorial aqui.
Coloque ambos no seu Diretório de Trabalho
Listas
Listas consistem em uma coleção de objetos, não necessariamente de mesma classe. Nelas podemos armazenar todos os outros objetos já vistos e recuperá-los pela indexação com [[. Como exemplo, vamos utilizar alguns objetos que já foram gerados.
minha_lista <- list(campo1 = campo1, media_alt = tapply(campo1$altura, campo1$idade, mean), matrix_ex = W)
str(minha_lista)## List of 3
## $ campo1 :'data.frame': 6 obs. of 6 variables:
## ..$ clone : chr [1:6] "GRA02" "URO01" "URO03" "GRA02" ...
## ..$ altura : num [1:6] 30.1 30.4 40 30.2 30.6 40.1
## ..$ diametro: num [1:6] 0.26 0.3 0.36 0.24 0.27 0.35
## ..$ idade : int [1:6] 3 4 5 3 4 5
## ..$ corte : logi [1:6] FALSE FALSE FALSE FALSE FALSE TRUE
## ..$ volume : num [1:6] 1.6 2.15 4.07 1.37 1.75 ...
## $ media_alt: num [1:3(1d)] 30.1 30.5 40
## ..- attr(*, "dimnames")=List of 1
## .. ..$ : chr [1:3] "3" "4" "5"
## $ matrix_ex: num [1:6, 1:2] 30.1 30.4 40 30.2 30.6 40.1 0.26 0.3 0.36 0.24 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:6] "GRA02" "URO01" "URO03" "GRA02" ...
## .. ..$ : chr [1:2] "altura" "diametro"Quero acessar o data.frame campo1
minha_lista[[1]] ## clone altura diametro idade corte volume
## 1 GRA02 30.1 0.26 3 FALSE 1.597287
## 2 URO01 30.4 0.30 4 FALSE 2.147760
## 3 URO03 40.0 0.36 5 FALSE 4.069440
## 4 GRA02 30.2 0.24 3 FALSE 1.365523
## 5 GRA01 30.6 0.27 4 FALSE 1.751131
## 6 URO01 40.1 0.35 5 TRUE 3.856116# ou
minha_lista$campo1## clone altura diametro idade corte volume
## 1 GRA02 30.1 0.26 3 FALSE 1.597287
## 2 URO01 30.4 0.30 4 FALSE 2.147760
## 3 URO03 40.0 0.36 5 FALSE 4.069440
## 4 GRA02 30.2 0.24 3 FALSE 1.365523
## 5 GRA01 30.6 0.27 4 FALSE 1.751131
## 6 URO01 40.1 0.35 5 TRUE 3.856116Para acessar uma coluna específica no data.frame campo1, que está dentro da minha_lista:
minha_lista[[1]][[3]]## [1] 0.26 0.30 0.36 0.24 0.27 0.35# ou
minha_lista[[1]]$diametro## [1] 0.26 0.30 0.36 0.24 0.27 0.35# ou
minha_lista$campo1$diametro## [1] 0.26 0.30 0.36 0.24 0.27 0.35Listas são muito úteis, por exemplo, quando vamos utilizar/gerar diversos objetos dentro de um loop.
Arrays
Este é um tipo de objeto que você provavelmente não irá utilizar agora no início, mas é bom saber da sua existência. São utilizados para armazenar dados com mais de duas dimensões. Por exemplo, se criarmos um array:
(meu_array <- array(1:24, dim = c(2,3,4)))## , , 1
##
## [,1] [,2] [,3]
## [1,] 1 3 5
## [2,] 2 4 6
##
## , , 2
##
## [,1] [,2] [,3]
## [1,] 7 9 11
## [2,] 8 10 12
##
## , , 3
##
## [,1] [,2] [,3]
## [1,] 13 15 17
## [2,] 14 16 18
##
## , , 4
##
## [,1] [,2] [,3]
## [1,] 19 21 23
## [2,] 20 22 24Teremos quatro matrizes com duas linhas e três colunas e os números de 1 a 24 estarão distribuídos nelas por colunas.
Se estiver adiantada/o em relação aos colegas, você já pode fazer os exercícios da Sessão 2, se não, faça-os em outro momento e nos envie dúvidas pelo fórum.
Importando e exportando dados
Os arquivos RData são exclusivos do R, clicando duas vezes no arquivo ou utilizando:
load("dia_1.RData")Você vai recuperar todos os objetos gerados até agora no tutorial. Para gerar esse arquivo RData, eu rodei todos os códigos daqui para cima e usei:
save.image(file = "dia_1.RData")Este comando salva tudo o que você tem no seu Ambiente Global (todos os objetos que aparecem ali no canto direito superior). Você pode também salvar somente um objeto com:
save(campo1, file = "campo1.RData")Se removermos ele do nosso Ambiente Global com:
rm(campo1) # Certifique-se que salvou o arquivo RData antes de removê-loPodemos facilmente obtê-lo novamente com:
load("campo1.RData")O formato RData é exclusivo para o R, ele é interessante de ser usado para fazer o que estamos fazendo, paramos a análise em um dia, vamos continuar em um outro e não queremos ter que rodar tudo de novo. Mas muitas vezes precisamos exportar nossos dados para outros programas, que exigem outros formatos, como, por exemplo, .txt ou .csv. Para isso utilizamos:
write.table(campo1, file = "campo1.txt", sep = ";", dec = ".", row.names = FALSE)
write.csv(campo1, file = "campo1.csv", row.names = TRUE)Obs: Você pode adquirir pacotes para exportar e importar dados com outros fomatos, como exemplo o pacote xlsx exporta e importa dados com formato do excel. O pacote vroom pode ser utilizado para compactar grandes tabelas.
Ao exportar, há diversas opções para a formatação do arquivo, é importante considerá-las se o arquivo for usado em outro sofware posteriormente.
Abra os arquivos gerados para visualizar sua formatação.
Esses arquivos podem ser lidos novamente pelo R, utilizando as funções e suas especificações:
campo1_txt <- read.table(file = "campo1.txt", sep=";", dec=".", header = TRUE)
campo1_csv <- read.csv(file = "campo1.csv")
head(campo1_txt)
head(campo1_csv)Agora que aprendemos a importar dados, vamos trabalhar com o conjunto gerado a partir dos participantes da disciplina Biometria de Marcadores - Turma de 2021.
A planilha com os dados está disponível no link abaixo, adicione-a ao seu diretório de trabalho ou indique o caminho da pasta ao importá-la para dentro do R, como a seguir.
Também podemos importar os dados diretamente do github, apontando o endereço web.
dados <- read.csv("https://github.com/Cristianetaniguti/Workshop_genetica_esalq/raw/master/SIFSC11/dados_alunos2021.csv")Aqui usaremos o argumento stringAsFactors que impede que o R transforme os vetores da tabela em fatores, os quais são mais difíceis de trabalhar. O argumento na.strings irá indicar como foram nomeados os dados perdidos.
dados <- read.csv(file = "dados_alunos2021.csv", stringsAsFactors = FALSE, na.strings="-", header = T, dec = ",")
head(dados)load("dados_alunos2021.RData")Vamos explorar a estrutura dos dados coletados:
str(dados)## 'data.frame': 116 obs. of 7 variables:
## $ Categoria. : chr "Doutorado" "Doutorado" "Mestrado" "Outro" ...
## $ Qual.seu.curso.de.graduação. : chr "Agronomia" "Agronomia" "Biotecnologia" "Licenciatura em Ciências Biológicas" ...
## $ Qual.seu.nível.de.conhecimento.sobre.Genética. : chr "Avançado" "Avançado" "Intermediário" "Intermediário" ...
## $ Qual.seu.nível.de.conhecimento.sobre.Estatística. : chr "Intermediário" "Intermediário" "Básico" "Intermediário" ...
## $ Qual.seu.nível.de.conhecimento.sobre..Genética.Estatística.. : chr "Intermediário" "Intermediário" "0" "Intermediário" ...
## $ Qual.a.latitude.do.local.onde.vai.estar.para.assistir.as.aulas..ATENÇÃO..use.notação.do.tipo..22.709556..com.decimais..e.confira.antes.de.postar..Não.é.usual.alunos.assistirem.as.aulas.no.meio.do.oceano..Use.o.Google.Maps.para.fazer.corretamente.: num -25.5 -22.7 -15.9 52.1 -22.8 ...
## $ Qual.a.longitude.do.local.onde.vai.estar.para.assistir.as.aulas..ATENÇÃO..use.notação.do.tipo..47.633866..com.decimais..e.confira.antes.de.postar..Não.é.usual.alunos.assistirem.as.aulas.no.meio.do.oceano..Use.o.Google.Maps.fazer.corretamente. : num -54.6 -47.6 -48.1 -106.7 -47.7 ...# também
dim(dados)## [1] 116 7Repare que nos nomes das colunas ainda estão as perguntas completas feitas no formulário, vamos alterar para nomes mais fáceis de trabalhar:
colnames(dados)## [1] "Categoria."
## [2] "Qual.seu.curso.de.graduação."
## [3] "Qual.seu.nível.de.conhecimento.sobre.Genética."
## [4] "Qual.seu.nível.de.conhecimento.sobre.Estatística."
## [5] "Qual.seu.nível.de.conhecimento.sobre..Genética.Estatística.."
## [6] "Qual.a.latitude.do.local.onde.vai.estar.para.assistir.as.aulas..ATENÇÃO..use.notação.do.tipo..22.709556..com.decimais..e.confira.antes.de.postar..Não.é.usual.alunos.assistirem.as.aulas.no.meio.do.oceano..Use.o.Google.Maps.para.fazer.corretamente."
## [7] "Qual.a.longitude.do.local.onde.vai.estar.para.assistir.as.aulas..ATENÇÃO..use.notação.do.tipo..47.633866..com.decimais..e.confira.antes.de.postar..Não.é.usual.alunos.assistirem.as.aulas.no.meio.do.oceano..Use.o.Google.Maps.fazer.corretamente."colnames(dados) = c("Ocupacao", "Graduacao", "Conhecimentos_Genetica", "Conhecimentos_Estatistica", "Conhecimento_Gen_Est",
"Latitude", "Longitude")
colnames(dados)## [1] "Ocupacao" "Graduacao"
## [3] "Conhecimentos_Genetica" "Conhecimentos_Estatistica"
## [5] "Conhecimento_Gen_Est" "Latitude"
## [7] "Longitude"str(dados)## 'data.frame': 116 obs. of 7 variables:
## $ Ocupacao : chr "Doutorado" "Doutorado" "Mestrado" "Outro" ...
## $ Graduacao : chr "Agronomia" "Agronomia" "Biotecnologia" "Licenciatura em Ciências Biológicas" ...
## $ Conhecimentos_Genetica : chr "Avançado" "Avançado" "Intermediário" "Intermediário" ...
## $ Conhecimentos_Estatistica: chr "Intermediário" "Intermediário" "Básico" "Intermediário" ...
## $ Conhecimento_Gen_Est : chr "Intermediário" "Intermediário" "0" "Intermediário" ...
## $ Latitude : num -25.5 -22.7 -15.9 52.1 -22.8 ...
## $ Longitude : num -54.6 -47.6 -48.1 -106.7 -47.7 ...Manipulando os dados
Agora usaremos os dados que temos para aprender diferentes comandos e funções do Ambiente R.
Primeiro, vamos verificar quantos alunos responderam as questões do formulario da disciplina, contando o número de linhas, para isso use a função nrow.
nrow(dados)## [1] 116Vamos então verificar se temos no nosso grupo pessoas que compartilham a mesma graduação e ocupação.
Podemos verificar isso facilmente com a função table, que indica a frequência de cada observação:
table(dados$Graduacao)
table(dados$Ocupacao)Estruturas condicionais
if e else
Para nossa próxima atividade com os dados, vamos primeiro entender como funcionam as estruturas if e else.
Nas funções condicionais if e else, estabelecemos uma condição para if, se ela for verdade a atividade será realizada, caso contrário (else) outra tarefa será. Como no exemplo:
if(2 >3){
cat("dois é maior que três")
} else {
cat("dois não é maior que três")
}## dois não é maior que três- Descubra a Nivel de Conhecimentos em Genetica (3ª coluna) da segunda pessoa que o respondeu (linha 2). Envie uma mensagem motivacional se ela possuir conhecimentos basicos, outra se a pessoa ja possuir conhecimento Intermediario ou Avançado (restante das respostas). (dica: o sinal
==se refere a “exatamente igual a”)
if(dados[2,3] =="Básico"){
cat("Força, Mendel acredita em você!")
} else {
cat("Mendel agradece a preferência")
}## Mendel agradece a preferênciaPodemos espeficiar mais do que uma condição repetindo a estrutura if else: Agora vamos estudar a ocupação das pessoas que responderam ao questionário.
if(dados[2,1] == "Mestrado"){
cat("Força, o quinto dia útil esta chegando!")
} else if (dados[2,1] == "Doutorado"){
cat("Assim como seus amigos do mestrado, acredite até sexta a bolsa cai!")
} else {
cat("Esse já tem a carteira assinada, que beleza!")
}## Assim como seus amigos do mestrado, acredite até sexta a bolsa cai!Mas repare que só é possível utilizar essas estruturas para um elemento individual do vetor ou em todo ele, se quisermos percorrer os elementos individualmente precisamos recorrer a outro recurso.
Estruturas de repetição
For
Esse recurso pode ser a função for, uma função muito utilizada e poderosa. Ela constitui uma estrutura de loop, pois irá aplicar a mesma atividade repetidamente até atingir uma determinada condição. Veja exemplos:
for(i in 1:10){
print(i)
}## [1] 1
## [1] 2
## [1] 3
## [1] 4
## [1] 5
## [1] 6
## [1] 7
## [1] 8
## [1] 9
## [1] 10test <- vector()
for(i in 1:10){
test[i] <- i+4
}
test## [1] 5 6 7 8 9 10 11 12 13 14Nos casos acima, i funciona como um índice que irá variar de 1 até 10 a operação determinada entre chaves.
Com essa estrutura, podemos repetir a operação realizada com as estruturas if e else para todo o vetor:
for(i in 1:nrow(dados)){
if(dados[i,1] == "Mestrado"){
print("Força, o quinto dia útil esta chegando!")
} else if (dados[i,1] == "Doutorado"){
print("Assim como seus amigos do mestrado, acredite até sexta a bolsa cai!")
} else {
print("Esse já tem a carteira assinada, que beleza!")
}
}Dica: Identação
Repare a diferença:
# Sem identação
for(i in 1:nrow(dados)){
if(dados[i,1] == "Mestrado"){
print("Força, o quinto dia útil esta chegando!")
} else if (dados[i,1] == "Doutorado"){
print("Assim como seus amigos do mestrado, acredite até sexta a bolsa cai!")
} else {
print("Esse já tem a carteira assinada, que beleza!")
}
}
# Com identação correta
for(i in 1:nrow(dados)){
if(dados[i,1] == "Mestrado"){
print("Força, o quinto dia útil esta chegando!")
} else if (dados[i,1] == "Doutorado"){
print("Assim como seus amigos do mestrado, acredite até sexta a bolsa cai!")
} else {
print("Esse já tem a carteira assinada, que beleza!")
}
}O editor de código do RStudio tem uma facilitação para identar códigos em R, selecione a área que deseja identar e aperte Ctrl+i.
Agora vamos trabalhar com a coluna 2, que possui a informação da graduação dos participantes. Repare que a função table retorna diferentes categorias que poderiam ser resumidas em apenas uma, como “Engenharia agronômica” e “Agronomia”. Vamos utilizar um loop para descobrir quais responderam áreas relacionadas à “Agro”. Depois, quais dessas respostas não foram “Agronomia” e pedir para que este participante modifique a resposta digitando apenas “Agronomia”.
Dica 1: Para identificar o padrão “Agro” podemos utilizar a função grepl. Dica 2: Podemos utilizar if dentro de if
# Exemplo do uso da função grepl
dados[,2]## [1] "Agronomia"
## [2] "Agronomia"
## [3] "Biotecnologia"
## [4] "Licenciatura em Ciências Biológicas"
## [5] "Agronomia"
## [6] "Agronomia"
## [7] "Zootecnia"
## [8] "Agronomia"
## [9] "Agronomia"
## [10] "Engenharia Florestal"
## [11] "Bacharel em Agronomia"
## [12] "Engenharia Agronômica"
## [13] "Engenharia Agrícola"
## [14] "Agronomia"
## [15] "Engenharia Agronômica"
## [16] "Ciências Biológicas"
## [17] "Agronomia"
## [18] "Engenheiro Florestal"
## [19] "Genética e Melhoramento de Plantas"
## [20] "Agronomia"
## [21] "Agronomia"
## [22] "Engenharia Agronômica"
## [23] "Agronomia"
## [24] "Agronomia"
## [25] "Ciências Biológicas"
## [26] "Agronomia"
## [27] "Ciências Biológicas"
## [28] "Engenharia Agronômica"
## [29] "Ciências Biológicas"
## [30] "Agronomia"
## [31] "Biotecnologia"
## [32] "Biologia"
## [33] "Agronomia"
## [34] "Biologia"
## [35] "Agronomia"
## [36] "Biologia"
## [37] "Agronomia"
## [38] "Melhoramento genético"
## [39] "FAEM - UFPEL"
## [40] "Agronomia"
## [41] "Ciências Biológicas"
## [42] "Ciências Biológicas"
## [43] "Agronomia"
## [44] "Biologia"
## [45] "Agronomia"
## [46] "Agronomia"
## [47] "Eng. Agronômica"
## [48] "Zootecnia"
## [49] "Agronomia"
## [50] "Biotecnologia"
## [51] "Agronomia"
## [52] "Agronomia"
## [53] "Agronomia"
## [54] "Agronomia"
## [55] "Agronomia"
## [56] "Agronomia"
## [57] "Agronomia"
## [58] "Genética e Melhoramento"
## [59] "Agronomia"
## [60] "Agronomia"
## [61] "Agronomia"
## [62] "Ciências bilógicas"
## [63] "Agronomia"
## [64] "Biologia"
## [65] "Genética e Melhoramento de Plantas"
## [66] "Agronomia"
## [67] "Agronomia"
## [68] "Agronomia"
## [69] "Agronomia"
## [70] "Agronomia"
## [71] "Agronomia"
## [72] "Agronomia"
## [73] "Agronomia"
## [74] "Engeheria em Biotecnologia Vegetal (Chile)"
## [75] "AGRONOMIA"
## [76] "Engenharia Agronômica"
## [77] "Agronomia"
## [78] "Agronomia"
## [79] "Biología"
## [80] "Agronomia"
## [81] "Ciências Biológicas"
## [82] "Agronomia"
## [83] "Agronomia"
## [84] "Agronomia"
## [85] "Agronomia"
## [86] "Agronomia"
## [87] "Ciências Biológicas"
## [88] "Agronomia"
## [89] "Agronomia"
## [90] "Agronomia"
## [91] "Agronomia"
## [92] "Agronomia"
## [93] "Engenharia Biotecnológica"
## [94] "Agronomia"
## [95] "Agronomia"
## [96] "Engenharia Florestal"
## [97] "Agronomia"
## [98] "Agronomia"
## [99] "Biologia"
## [100] "Biologia"
## [101] "Ciências Biológicas"
## [102] "Agronomia"
## [103] "Agronomia"
## [104] "Agronomia"
## [105] "Engenheiro Agrônomo"
## [106] "Genética e Melhoramento"
## [107] "Agronomia"
## [108] "Agronomia"
## [109] "Biotecnologia"
## [110] "Agronomia"
## [111] "Agronomia"
## [112] "Biotecnologia"
## [113] "Engenharia Agronômica / Licenciatura em Ciências Agrárias"
## [114] "Ciencias Biologicas"
## [115] "Agronomia"
## [116] "Eng. Florestal"grepl("Agro", dados[,2]) # Quais linhas contém os caracteres "Agro"## [1] TRUE TRUE FALSE FALSE TRUE TRUE FALSE TRUE TRUE FALSE TRUE TRUE
## [13] FALSE TRUE TRUE FALSE TRUE FALSE FALSE TRUE TRUE TRUE TRUE TRUE
## [25] FALSE TRUE FALSE TRUE FALSE TRUE FALSE FALSE TRUE FALSE TRUE FALSE
## [37] TRUE FALSE FALSE TRUE FALSE FALSE TRUE FALSE TRUE TRUE TRUE FALSE
## [49] TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE
## [61] TRUE FALSE TRUE FALSE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
## [73] TRUE FALSE FALSE TRUE TRUE TRUE FALSE TRUE FALSE TRUE TRUE TRUE
## [85] TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE FALSE
## [97] TRUE TRUE FALSE FALSE FALSE TRUE TRUE TRUE FALSE FALSE TRUE TRUE
## [109] FALSE TRUE TRUE FALSE TRUE FALSE TRUE FALSEdados[grepl("Agro", dados[,2]),2]## [1] "Agronomia"
## [2] "Agronomia"
## [3] "Agronomia"
## [4] "Agronomia"
## [5] "Agronomia"
## [6] "Agronomia"
## [7] "Bacharel em Agronomia"
## [8] "Engenharia Agronômica"
## [9] "Agronomia"
## [10] "Engenharia Agronômica"
## [11] "Agronomia"
## [12] "Agronomia"
## [13] "Agronomia"
## [14] "Engenharia Agronômica"
## [15] "Agronomia"
## [16] "Agronomia"
## [17] "Agronomia"
## [18] "Engenharia Agronômica"
## [19] "Agronomia"
## [20] "Agronomia"
## [21] "Agronomia"
## [22] "Agronomia"
## [23] "Agronomia"
## [24] "Agronomia"
## [25] "Agronomia"
## [26] "Agronomia"
## [27] "Eng. Agronômica"
## [28] "Agronomia"
## [29] "Agronomia"
## [30] "Agronomia"
## [31] "Agronomia"
## [32] "Agronomia"
## [33] "Agronomia"
## [34] "Agronomia"
## [35] "Agronomia"
## [36] "Agronomia"
## [37] "Agronomia"
## [38] "Agronomia"
## [39] "Agronomia"
## [40] "Agronomia"
## [41] "Agronomia"
## [42] "Agronomia"
## [43] "Agronomia"
## [44] "Agronomia"
## [45] "Agronomia"
## [46] "Agronomia"
## [47] "Agronomia"
## [48] "Engenharia Agronômica"
## [49] "Agronomia"
## [50] "Agronomia"
## [51] "Agronomia"
## [52] "Agronomia"
## [53] "Agronomia"
## [54] "Agronomia"
## [55] "Agronomia"
## [56] "Agronomia"
## [57] "Agronomia"
## [58] "Agronomia"
## [59] "Agronomia"
## [60] "Agronomia"
## [61] "Agronomia"
## [62] "Agronomia"
## [63] "Agronomia"
## [64] "Agronomia"
## [65] "Agronomia"
## [66] "Agronomia"
## [67] "Agronomia"
## [68] "Agronomia"
## [69] "Agronomia"
## [70] "Agronomia"
## [71] "Agronomia"
## [72] "Agronomia"
## [73] "Engenharia Agronômica / Licenciatura em Ciências Agrárias"
## [74] "Agronomia"for(i in 1:nrow(dados)){
if(grepl("Agro", dados[i,2])){
if(dados[i,2] != "Agronomia"){
print("Por favor, substitua sua resposta por Agronomia.")
}
}
}## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."Para que seja possível imprimir conteúdo de objetos durante o loop, usamos a função cat, ela não separa cada resposta em uma linha, precisamos colocar o \n indicando a quebra de linha.
Agora vamos nos mesmos homogeneizar essas informações. Podemos armazenar em uma variável a posição das linhas identificadas, então corrigiremos manualmente somente essas:
homog <- vector()
for(i in 1:nrow(dados)){
if(grepl("Agro", dados[i,2])){
if(dados[i,2] != "Agronomia"){
print("Por favor, substitua sua resposta por Agronomia.")
homog <- c(homog, i)
}
}
}## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."
## [1] "Por favor, substitua sua resposta por Agronomia."homog## [1] 11 12 15 22 28 47 76 113Como você faria para corrigir esses elementos errados? Tente!
‘dados[homog, 2] <- “Agronomia”’
While
Nesse tipo de estrutura de repetição a tarefa será realizada até que seja atingida determinada condição.
x <- 1
while(x < 5){
x <- x + 1
cat(x)
}## 2345É muito importante que nessa estrutura a condição seja atingida, caso contrário o loop irá funcionar infinitamente e você terá que interrompê-lo por meios externos, como, se este utilizando RStudio, clicar no simbolo em vermelho no canto direito superior da janela do console, ou apertar Ctrl+C no console.
Não é muito difícil disso acontecer, basta um pequeno erro como:
x <- 1
while(x < 5){
x + 1
cat(x)
}Aqui podemos utilizar os comandos break e next para atender a outras condições, como:
x <- 1
while(x < 5){
x <- x + 1
if(x==4) break
cat(x)
}## 23x <- 1
while(x < 5){
x <- x + 1
if(x==4) next
cat(x)
}## 235Repeat
Esta estrutura também exige uma condição de parada, mas esta condição é necessariamente colocada dentro do bloco de código com o uso do break. Ela então repete o bloco de código até a condição o interrompa.
x <- 1
repeat{
x <- x+1
cat(x)
if(x==4) break
}## 234Loops dentro de loops
É possível também utilizarmos estruturas de repetição dentro de estruturas de repetição. Por exemplo, se quisermos trabalhar tanto nas colunas como nas linhas de uma matrix.
# Criando uma matrix vazia
ex_mat <- matrix(nrow=10, ncol=10)
# cada número dentro da matrix será o produto no índice da coluna pelo índice da linha
for(i in 1:dim(ex_mat)[1]) {
for(j in 1:dim(ex_mat)[2]) {
ex_mat[i,j] = i*j
}
}Outro exemplo de uso:
var1 <- c("fertilizante1", "fertilizante2")
var2 <- c("ESS", "URO", "GRA")
w <- 1
for(i in var1){
for(j in var2){
nome_arquivo <- paste0(i,"_planta_",j,".txt")
arquivo <- data.frame("bloco" = "fake_data", "tratamento" ="fake_data")
write.table(arquivo, file = nome_arquivo)
w <- w + 1
}
}
# Verifique seu diretorio de trabalho, arquivos devem ter sido geradosFizemos um vídeo com mais detalhes sobre loops no R (quando cursamos a disciplina de Biometria de Marcadores), aumentem nossa quantidade de views e likes por lá.
Se estiver adiantada/o em relação aos colegas, você já pode fazer os exercícios da Sessão 3, se não, faça-os em outro momento e nos envie dúvidas pelo fórum.
Algumas dicas:
- Cuidado ao rodar o mesmo comando mais de uma vez, algumas variáveis podem não ser mais como eram antes. Para que o comando funcione da mesma forma é necessário que os objetos de entrada estejam da forma como você espera.
- Lembrem-se que
=é para definir objetos e==é o sinal de igualdade. - Nas estruturas condicionais e de repetição, lembrem-se que é necessário manter a sintaxe esperada: If(){} e for(i in 1:10){}. No for, podemos trocar a letra que será o índice, mas é sempre necessário fornecer uma sequência de inteiros ou caracteres.
- Usar identação ajuda a visualizar o começo e fim de cada estrutura de código e facilita o abrir e fechar de chaves. Identação são aqueles espaços que usamos antes da linha, como:
# Criando uma matrix vazia
ex_mat <- matrix(nrow=10, ncol=10)
# cada número dentro da matrix será o produto no índice da coluna pelo índice da linha
for(i in 1:dim(ex_mat)[1]) { # Primeiro nível, não tem espaço
for(j in 1:dim(ex_mat)[2]) { # Segundo nível tem um espaço (tab)
ex_mat[i,j] = i*j # Terceiro nível tem dois espaços
} # Fechei o segundo nível
} # Fechei o primeiro nívelElaboração de gráficos simples
Para outros dados coletados, vamos gerar alguns gráficos simples utilizando as funções básicas do R. Existem pacotes como o ggplot2, plotly e shiny que possuem ferramentas muito poderosas para construção de gráficos, mas exigem um pouco mais de tempo para aprendizagem de sua sintaxe.
Os tipos mais comuns já possuem funções próprias, mas outros gráficos podem ser customizados de acordo com a necessidade do usuário. Vamos iniciar com um simples gráfico de frequências (ou histograma) para os dados de latitude.
hist(dados$Latitude)
Vamos adicionar alguns argumentos para dar uma personalizada:
breakspara definir os intervalos para cada barra;
Agora tente fazer o histograma para a Longitude, aproveite para tentar alterar alguns parâmetros. Em seguida, serão apresentados outros gráficos que poderão ser utilizados.
Salvar gráficos
Os gráficos podem ser salvos através dos menus disponíveis no RStudio, ou através de funções que permitem salvar em formatos específicos. Algumas delas são: pdf(); png(); jpeg(); bitmap(). De maneira geral, o parâmetro primordial é fornecer o nome do arquivo que será gerado (contendo sua extensão). Após abrir a função gráfica, deve-se gerar o gráfico de interesse. Por fim, utiliza-se o comando dev.off() para que saída gráfica volte para o console.
png(filename = "hist_rbase.png")
hist(dados$Latitude)
dev.off()
png(filename = "hist_rbase.png", width = 1500, height = 1500, res= 300)
hist(dados$Latitude)
dev.off()Agora, gere um gráfico e salve-o no formato de seu interesse. Em seguida, crie diversos gráficos dentro de uma mesma função gráfica e estude a saída.
Um exemplo didático de aplicação em genética
Aqui vamos exercitar um pouco do que aprendemos de R e já introduzir uma das muitas aplicações em genética.
Para melhor didática, este é um exemplo hipotético, a partir de dados simulados para uma população de retrocruzamento.
Mas vamos imaginar que se trata de um experimento de avaliação da resistência de diferentes cultivares de batata contra um fungo patógeno.
Aqui temos a média ajustadas das medida do tamanho da lesão nas folhas de cada cultivar inoculadas.
fenotipos <- read.csv("https://raw.githubusercontent.com/Cristianetaniguti/Workshop_genetica_esalq/gh-pages/4biotec/fenotipos.csv",
stringsAsFactors = F)
str(fenotipos)## 'data.frame': 20 obs. of 3 variables:
## $ X : int 1 2 3 4 5 6 7 8 9 10 ...
## $ cultivar : chr "Cult1" "Cult2" "Cult3" "Cult4" ...
## $ tamanho.da.lesÃ.o: num 6.55 3.94 2.79 3.68 1.57 ...Ajuste os dados deixando apenas um coluna referente ao tamanho da lesão e o nome das cultivares como nome das linhas.
rownames(fenotipos) <- fenotipos$cultivar
fenotipos <- fenotipos[,-c(1,2)]E aqui temos os genótipos de cada uma das cultivares. Podemos obtê-los via marcadores moleculares, como os SNPs.
Aa = 1
aa = 0genotipos <- read.csv("https://raw.githubusercontent.com/Cristianetaniguti/Workshop_genetica_esalq/gh-pages/4biotec/genotipos.csv",
stringsAsFactors = F)
head(genotipos)## X MK1 MK2 MK3 MK4 MK5 MK6 MK7 MK8 MK9 MK10 MK11 MK12 MK13 MK14 MK15 MK16
## 1 Cult1 1 0 0 0 1 1 0 1 1 1 1 1 1 0 1 1
## 2 Cult2 1 1 1 0 1 1 0 0 1 0 0 0 0 1 1 1
## 3 Cult3 1 1 0 0 1 1 0 0 0 1 0 1 0 0 1 1
## 4 Cult4 1 1 0 0 0 0 0 1 0 1 0 0 0 0 1 1
## 5 Cult5 0 1 1 0 0 0 0 0 0 0 1 1 0 0 0 1
## 6 Cult6 1 0 1 0 1 0 0 0 1 1 0 0 0 1 1 0
## MK17 MK18 MK19 MK20 MK21 MK22 MK23 MK24 MK25
## 1 0 1 1 1 1 1 1 1 1
## 2 1 0 0 1 1 1 1 1 1
## 3 1 0 0 0 0 0 1 1 1
## 4 1 0 0 1 1 0 1 1 0
## 5 0 0 0 0 0 0 1 1 0
## 6 1 0 0 0 0 0 1 1 1Observe que a primeira coluna da matrix é o nome da cultivar. Vamos fazer com que ela vire o nome das linhas.
rownames(genotipos) <- genotipos[,1]
genotipos <- genotipos[,-1]
head(genotipos)## MK1 MK2 MK3 MK4 MK5 MK6 MK7 MK8 MK9 MK10 MK11 MK12 MK13 MK14 MK15 MK16
## Cult1 1 0 0 0 1 1 0 1 1 1 1 1 1 0 1 1
## Cult2 1 1 1 0 1 1 0 0 1 0 0 0 0 1 1 1
## Cult3 1 1 0 0 1 1 0 0 0 1 0 1 0 0 1 1
## Cult4 1 1 0 0 0 0 0 1 0 1 0 0 0 0 1 1
## Cult5 0 1 1 0 0 0 0 0 0 0 1 1 0 0 0 1
## Cult6 1 0 1 0 1 0 0 0 1 1 0 0 0 1 1 0
## MK17 MK18 MK19 MK20 MK21 MK22 MK23 MK24 MK25
## Cult1 0 1 1 1 1 1 1 1 1
## Cult2 1 0 0 1 1 1 1 1 1
## Cult3 1 0 0 0 0 0 1 1 1
## Cult4 1 0 0 1 1 0 1 1 0
## Cult5 0 0 0 0 0 0 1 1 0
## Cult6 1 0 0 0 0 0 1 1 1Como já vimos, gráficos são ferramentas interessantes para melhor visualização dos dados. Aqui vamos construir um heatmap, cada genótipo receberá uma cor diferente. Assim poderemos visualizar melhor a distribuição deles na população estudada.
heatmap(genotipos, Colv = NA, Rowv = NA)Opa! Aconteceu algo errado, né? É importante prestar atenção no que o erro diz, essas dicas são valiosas. O problema aqui é que a função só recebe objetos matriz e fornecemos para ela um data.frame.
genotipos <- as.matrix(genotipos)
heatmap(genotipos, Colv = NA, Rowv = NA)
Para que o gráfico fique intuitivo, precisamos fazer algumas alterações. Repare que embora a gente tenha somente genótipos, existem mais cores no gráfico. Vamos corrigir isso.
geno.cores <- heat.colors(2)
heatmap(genotipos, Colv = NA, Rowv = NA, col=geno.cores)
head(genotipos)## MK1 MK2 MK3 MK4 MK5 MK6 MK7 MK8 MK9 MK10 MK11 MK12 MK13 MK14 MK15 MK16
## Cult1 1 0 0 0 1 1 0 1 1 1 1 1 1 0 1 1
## Cult2 1 1 1 0 1 1 0 0 1 0 0 0 0 1 1 1
## Cult3 1 1 0 0 1 1 0 0 0 1 0 1 0 0 1 1
## Cult4 1 1 0 0 0 0 0 1 0 1 0 0 0 0 1 1
## Cult5 0 1 1 0 0 0 0 0 0 0 1 1 0 0 0 1
## Cult6 1 0 1 0 1 0 0 0 1 1 0 0 0 1 1 0
## MK17 MK18 MK19 MK20 MK21 MK22 MK23 MK24 MK25
## Cult1 0 1 1 1 1 1 1 1 1
## Cult2 1 0 0 1 1 1 1 1 1
## Cult3 1 0 0 0 0 0 1 1 1
## Cult4 1 0 0 1 1 0 1 1 0
## Cult5 0 0 0 0 0 0 1 1 0
## Cult6 1 0 0 0 0 0 1 1 1Agora eu quero visualizar a relação desses genótipos com os fenótipos. Para isso, vou acrescentar uma coluna ao lado com cores conforme os valores dos fenótipos. Ou seja, farei outro gradiente com cores frias para valores baixos do tamanho da lesão e cores fortes para valores altos.
fenotipos_ordenado <- sort(fenotipos)
desordem <- match(fenotipos, fenotipos_ordenado)
fenotipos## [1] 6.5524741 3.9354826 2.7900899 3.6761125 1.5747197 2.0330726
## [7] 8.0846823 -0.8723574 8.0466798 3.8435632 5.4567962 2.2027866
## [13] 2.9866094 3.6863699 6.5323953 4.7906945 5.2409508 5.7539092
## [19] 3.3656062 6.3934879fenotipos_ordenado## [1] -0.8723574 1.5747197 2.0330726 2.2027866 2.7900899 2.9866094
## [7] 3.3656062 3.6761125 3.6863699 3.8435632 3.9354826 4.7906945
## [13] 5.2409508 5.4567962 5.7539092 6.3934879 6.5323953 6.5524741
## [19] 8.0466798 8.0846823desordem## [1] 18 11 5 8 2 3 20 1 19 10 14 4 6 9 17 12 13 15 7 16colfunc <- colorRampPalette(c("royalblue", "red"))
cores <- colfunc(20)[desordem]
heatmap(genotipos, Colv = NA, Rowv = NA, col=heat.colors(2), RowSideColors = cores)
Uma bagunça! Não consigo tirar nenhuma conclusão. Mas.. e se eu ordenar conforme o valor do fenótipos?
ordem <- match(fenotipos_ordenado, fenotipos)
heatmap(genotipos[ordem,], Rowv = NA, Colv = NA, col=geno.cores, RowSideColors = cores[ordem])
Quase conseguimos ver uma tendência.
Como são dados simulados, sabemos que os marcadores de 6 a 10 e de 18 a 20 estão relacionados com a característica. Vamos veer o que acontece se separamos eles.
heatmap(genotipos[ordem,c(6:10,18:20)], Rowv = NA, Colv = NA, col=geno.cores, RowSideColors = cores[ordem])
Não é interessante?!? E ai? Qual genótipo é vantajoso ter nesses locos para o plantio de batata?
E os marcadores não ligados?
heatmap(genotipos[ordem,-c(6:10,17:20)], Rowv = NA, Colv = NA, col=geno.cores, RowSideColors = cores[ordem])
A bagunça continua! Indicando uma aleatoriedade, ou seja, não existe uma correlação entre o genótipo e o fenótipo.
Mas uma análise visual apenas instiga nossa intuição, precisamos de uma análise estatística para realmente verificarmos. Existem muitos métodos estatísticos desenvolvidos para isso, muito robustos que levam em conta vários fatores importante como o desequilibrio de ligação entre os marcadores e a influência da existência de efeitos próximos. Mas por agora faremos a mais simples das análises, já muito obsoleta para esse propósito, mas servirá para o nosso aprendizado. Faremos uma análise de variância para cada marcador individualmente.
mk <- genotipos[,11]
mk <- gsub(pattern = 0, replacement = "aa", x = mk)
mk <- gsub(pattern = 1, replacement = "Aa", x = mk)
mk <- gsub(pattern = 2, replacement = "AA", x = mk)
mk <- as.factor(mk)
modelo <- lm(fenotipos ~ mk)
summary(modelo)##
## Call:
## lm(formula = fenotipos ~ mk)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.7228 -1.1207 -0.1692 1.3946 4.2343
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.8504 0.6183 6.228 7.09e-06 ***
## mkAa 1.2951 1.0451 1.239 0.231
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.229 on 18 degrees of freedom
## Multiple R-squared: 0.07862, Adjusted R-squared: 0.02743
## F-statistic: 1.536 on 1 and 18 DF, p-value: 0.2312mk <- genotipos[,10]
mk <- gsub(pattern = 0, replacement = "aa", x = mk)
mk <- gsub(pattern = 1, replacement = "Aa", x = mk)
mk <- gsub(pattern = 2, replacement = "AA", x = mk)
mk <- as.factor(mk)
modelo <- lm(fenotipos ~ mk)
summary(modelo)##
## Call:
## lm(formula = fenotipos ~ mk)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.7812 -0.8631 0.3829 1.0418 2.6397
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.9088 0.6346 4.584 0.00023 ***
## mkAa 2.5361 0.8557 2.964 0.00832 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.904 on 18 degrees of freedom
## Multiple R-squared: 0.3279, Adjusted R-squared: 0.2906
## F-statistic: 8.784 on 1 and 18 DF, p-value: 0.008315Não queremos ficar digitando este código 25 vezes, né? Vamos aplicar nossos conhecimentos de estruturas de repetição (loops)
p <- vector()
for(i in 1:dim(genotipos)[2]){
mk <- genotipos[,i]
mk <- gsub(pattern = 0, replacement = "aa", x = mk)
mk <- gsub(pattern = 1, replacement = "Aa", x = mk)
mk <- gsub(pattern = 2, replacement = "AA", x = mk)
mk <- as.factor(mk)
modelo <- lm(fenotipos ~ mk)
f <- summary(modelo)$fstatistic
p[i] <- pf(f[1],f[2],f[3],lower.tail=F)
}
which(p < 0.05)## [1] 6 7 8 9 10 18 20 21Aplicações de pacotes
Existem diversos pacotes disponíveis para variadas aplicações. Utilizaremos o ggplot2, que está disponível no repositório oficial do R, o CRAN. Portanto para instalá-lo:
install.packages("ggplot2")Depois disso é necessário recrutá-lo toda vez que iniciar um sessão no R:
library("ggplot2")O ggplot2 é um pacote que permite a construção de gráficos estatísticos, suas funcionalidades vão muito além do que está disponível nos gráficos básicos do R. Saiba mais sobre ele neste link. Encontre aqui vários exemplos de como ele pode ser usado.
Pratique gerando relatórios no RStudio
Utilize o R no seu dia-a-dia para ir praticando a linguagem. Além das recomendações contidas na primeira apresentação, recomendamos também dar uma olhada em como gerar documentos em pdf e html usando a Markdown. Utilizamos essa metodologia para gerar este tutorial e outras apresentações do curso. Pode ser muito prático no dia-a-dia!
Para utilizar, será necessário a instalação de outros pacotes. Um deles é o próprio rmarkdown:
install.packages("rmarkdown")library(rmarkdown)Agora crie um arquivo .Rmd utilizando as facilidades do RStudio, clique no ícone com símbolo + no canto superior esquerdo. Escolha o opção R Markdown. Dê um título ao seu arquivo e escolha a opção html. Ao fazer isso, o RStudio já coloca um template inicial, ja com um cabeçalho:
---
title: "Teste"
author: "Eu"
date: "June 5, 2018"
output: html_document
---
Este é o mais simples possível, você pode otimizá-lo de diversas maneiras. Saiba mais aqui.
O template inicial também traz alguns exemplos de sintaxe do markdown. Observe que utilizando # para títulos de sessões, ## para um nível inferior (subtitulos) e assim por diante. Palavras em negrito são escritas em meia a dois * e existem diversas outras especificações para essa sintaxe. Veja mais sobre ela aqui.
Para compilar o código, basta clicar em Knit. Ele irá pedir para que o arquivo .Rmd seja salvo com algum nome em algum lugar.
O markdown também é capaz de entender diretamente a linguagem html, também a css e latex. Para essa última, o latex precisa estar instalado e todas suas dependências.
Existem alguns pacotes que fornecem templates mais robustos para produção de htmls. Para esse tutorial utilizando o pacote rmdformats e personalizamos suas cores. Experimente:
install.packages("rmdformats")Agora faça o mesmo procedimento, clique no +, escolha R Markdown e, antes de escolher um título, mude para From Template, escolha o HTML readthedown template. Copie e cole o seguinte texto e aperte Knit.
# Teste1
Isso aqui é um teste só para dar uma olhada no template
## Testinho
Subsessão
* Item
**negrito**
*itálico*
fiz um [link](https://GENt-esalq.github.io/)!
Saiba mais no tutorial do R-bloggers, que acreditamos ser um bom começo! Acesse aqui.
Extra
Algumas ferramentas básicas de análise de dados
Claramente a análise de dados é algo muito específico de cada conjunto de dados e interesses, aqui apenas exemplificaremos o uso de algumas funções do R em duas situações específicas.
Avaliando o clima de Londrina
Vamos utilizar outro conjunto de dados para realizarmos mais avaliações utilizando a função lm. Acesse o conjunto clima_lond:
load("clima_lond.RData")Para obter os dados de precipitação da cidade de Londrina no primeiro semestre de 2017. Vamos utilizar as funções tapply e lm para avaliar os dados.
Algumas avaliações descritivas podem ser feitas pelo uso do tapply e de gráficos. A função summary também fornece informações gerais do conjunto. É possível usá-la em conjunto com o tapply.
# Verificando se as variáveis categórias estão como fatores
str(clima_lond)## 'data.frame': 181 obs. of 4 variables:
## $ dia : Factor w/ 31 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ Mes : Factor w/ 6 levels "Abril","Fevereiro",..: 3 3 3 3 3 3 3 3 3 3 ...
## $ prec.mm: Factor w/ 56 levels "0","0.2","0.3",..: 55 52 26 5 21 1 4 1 1 4 ...
## $ Data : Factor w/ 181 levels "2017-01-01","2017-01-02",..: 1 2 3 4 5 6 7 8 9 10 ...clima_lond$dia <- as.factor(clima_lond$dia)
# A precipitação nesse caso é uma variável contínua, nao categórica, para transformá-la use:
clima_lond$prec.mm <- as.numeric(as.character(clima_lond$prec.mm))
# Já com o tapply podemos ver as diferenças
tapply(clima_lond$prec.mm, clima_lond$Mes, summary)## $Abril
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 4.653 0.650 47.400
##
## $Fevereiro
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 0.00 0.00 3.65 2.85 34.80
##
## $Janeiro
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.600 8.532 6.600 59.600
##
## $Junho
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 3.573 0.500 30.600
##
## $Maio
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 7.406 7.000 48.200
##
## $`Março`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 3.223 1.200 26.200Repare que os níveis aparecem em ordem alfabética e não conforme o tempo, alteramos isso com:
levels(clima_lond$Mes)## [1] "Abril" "Fevereiro" "Janeiro" "Junho" "Maio" "Março"# A função match vai indicar a posição dos elementos do primeiro vetor no segundo vetor
pos <- match(c("Janeiro", "Fevereiro", "Março", "Abril", "Maio", "Junho"), levels(clima_lond$Mes))
pos## [1] 3 2 NA 1 5 4# Aqui utilizamos as posições obtidas para reordenar os meses nos níveis do fator
clima_lond$Mes = factor(clima_lond$Mes,
levels(clima_lond$Mes)[pos])
# Refazendo
tapply(clima_lond$prec.mm, clima_lond$Mes, summary)## $Janeiro
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.600 8.532 6.600 59.600
##
## $Fevereiro
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 0.00 0.00 3.65 2.85 34.80
##
## $Abril
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 4.653 0.650 47.400
##
## $Maio
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 7.406 7.000 48.200
##
## $Junho
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.000 0.000 0.000 3.573 0.500 30.600Podemos também avaliar visualmente através de gráficos:
ggplot(clima_lond) +
geom_point(aes(x = Mes, y = prec.mm)) +
labs(title = "Precipitação x Mês", x = "Meses do ano de 2017", y = "Precipitação em mm")
barplot(tapply(clima_lond$prec.mm, clima_lond$Mes, sum),
main="Total Mensal",
xlab = "Meses do ano de 2017",
ylab = "Precipitação em mm")
# Vamos fazer um gráfico de barras que mostre a soma de precipitação em cada mês
stat_lond <- tapply(clima_lond$prec.mm, clima_lond$Mes, sum)
str(stat_lond)## num [1:5(1d)] 264 102 140 230 107
## - attr(*, "dimnames")=List of 1
## ..$ : chr [1:5] "Janeiro" "Fevereiro" "Abril" "Maio" ...# Repare que o resultado do tapply é um vetor com seus elementos nomeados com cada mês
# Para usar o ggplot nesse caso, precisamos elaborar um data.frame cujas colunas sejam as variáveis que utilizaremos
stat_lond_edit <- data.frame("mes" = factor(names(stat_lond),
levels = c("Janeiro", "Fevereiro", "Março", "Abril", "Maio", "Junho")), "soma" = stat_lond)
str(stat_lond_edit)## 'data.frame': 5 obs. of 2 variables:
## $ mes : Factor w/ 6 levels "Janeiro","Fevereiro",..: 1 2 4 5 6
## $ soma: num 264 102 140 230 107ggplot(stat_lond_edit) +
geom_bar(aes(x = mes, y = soma), stat="identity") +
labs(title= "Soma da precipitação por mês", x = "Meses do ano de 2017", y = "Soma da precipitação")
Vamos então realizar uma análise de variância a fim de verificar se há diferenças significativas entre os meses. As funções lm e summary diferem apenas na forma de apresentar os resultados. O p-valor nos indica se podemos considerar diferenças do peso conforme o gênero.
mod <- lm(prec.mm ~ Mes, data = clima_lond)
summary(mod)##
## Call:
## lm(formula = prec.mm ~ Mes, data = clima_lond)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.532 -7.406 -3.650 -1.125 51.068
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.532 2.181 3.912 0.00014 ***
## MesFevereiro -4.882 3.166 -1.542 0.12524
## MesAbril -3.879 3.110 -1.247 0.21434
## MesMaio -1.126 3.085 -0.365 0.71565
## MesJunho -4.959 3.110 -1.594 0.11301
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12.14 on 145 degrees of freedom
## (31 observations deleted due to missingness)
## Multiple R-squared: 0.02836, Adjusted R-squared: 0.00156
## F-statistic: 1.058 on 4 and 145 DF, p-value: 0.3795modaov <- aov(prec.mm ~ Mes, data = clima_lond)
summary(modaov)## Df Sum Sq Mean Sq F value Pr(>F)
## Mes 4 624 156.1 1.058 0.38
## Residuals 145 21384 147.5
## 31 observations deleted due to missingnessPodemos também fazer um teste de médias para diferenciar a precipitação no decorrer dos meses. Aqui utilizaremos o método de Tukey:
tukey.test <- TukeyHSD(x=modaov, 'Mes', conf.level=0.95)
tukey.test## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = prec.mm ~ Mes, data = clima_lond)
##
## $Mes
## diff lwr upr p adj
## Fevereiro-Janeiro -4.88225806 -13.628251 3.863735 0.5369940
## Abril-Janeiro -3.87892473 -12.470374 4.712524 0.7236982
## Maio-Janeiro -1.12580645 -9.646543 7.394930 0.9961749
## Junho-Janeiro -4.95892473 -13.550374 3.632524 0.5033691
## Abril-Fevereiro 1.00333333 -7.811566 9.818232 0.9978586
## Maio-Fevereiro 3.75645161 -4.989541 12.502445 0.7591298
## Junho-Fevereiro -0.07666667 -8.891566 8.738232 0.9999999
## Maio-Abril 2.75311828 -5.838331 11.344567 0.9019673
## Junho-Abril -1.08000000 -9.741585 7.581585 0.9969451
## Junho-Maio -3.83311828 -12.424567 4.758331 0.7324452Avaliando experimento de café
Agora, trabalharemos com dados de um experimento de café. Acesse aqui:
- Arquivo cafe.txt
O experimento trata-se de dados em blocos completos casualizados de 10 progênies de café. Nele a coluna rep refere-se à repetição, prog indica o indivíduo da progênie (prog) e colheita indica a colheita.
data <- read.table("cafe.txt", h = TRUE, sep = "\t", dec = ",")
str(data)## 'data.frame': 120 obs. of 4 variables:
## $ rep : int 1 1 1 1 1 1 1 1 1 1 ...
## $ prog : int 1 2 3 4 5 6 7 8 9 10 ...
## $ colheita: int 1 1 1 1 1 1 1 1 1 1 ...
## $ prod : num 4.3 13.2 1.59 2.76 11.38 ...Não esqueça que é necessário que o arquivo esteja no seu ambiente de trabalho ou que você especifique o caminho completo para que o R o encontre!
Para essa análise de dados, nossa variável resposta é a produção (prod) e a repetição (rep), a progênie (prog) e a colheita (colheita) serão fatores no nosso modelo, identificados por seus níveis.
# Transformar em fator
data$rep <- as.factor(data$rep)
data$prog <- as.factor(data$prog)
data$colheita <- as.factor(data$colheita)
str(data)## 'data.frame': 120 obs. of 4 variables:
## $ rep : Factor w/ 4 levels "1","2","3","4": 1 1 1 1 1 1 1 1 1 1 ...
## $ prog : Factor w/ 10 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ colheita: Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ prod : num 4.3 13.2 1.59 2.76 11.38 ...# Outra opção
data <- transform(data, rep = factor(rep), prog = factor(prog), colheita = factor(colheita))
str(data)## 'data.frame': 120 obs. of 4 variables:
## $ rep : Factor w/ 4 levels "1","2","3","4": 1 1 1 1 1 1 1 1 1 1 ...
## $ prog : Factor w/ 10 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ colheita: Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ prod : num 4.3 13.2 1.59 2.76 11.38 ...Vamos analisar somente os dados referentes à primeira colheita. Para isso, podemos obter o subconjunto referente a ela:
# Indexar primeita colheita
Colheita_1 <- subset(data, colheita == 1)
str(Colheita_1)## 'data.frame': 40 obs. of 4 variables:
## $ rep : Factor w/ 4 levels "1","2","3","4": 1 1 1 1 1 1 1 1 1 1 ...
## $ prog : Factor w/ 10 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ colheita: Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ prod : num 4.3 13.2 1.59 2.76 11.38 ...Repare que, ao fazer o subconjunto, os três níveis do fator colheita ainda são mantidos, embora agora tenhamos apenas um. Isso pode ser um problema para a nossa análise. Portanto, devemos remover os níveis excedentes:
# Droplevels
Colheita_1 <- droplevels(subset(data, colheita == 1))
str(Colheita_1)## 'data.frame': 40 obs. of 4 variables:
## $ rep : Factor w/ 4 levels "1","2","3","4": 1 1 1 1 1 1 1 1 1 1 ...
## $ prog : Factor w/ 10 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ colheita: Factor w/ 1 level "1": 1 1 1 1 1 1 1 1 1 1 ...
## $ prod : num 4.3 13.2 1.59 2.76 11.38 ...Em seguida, podemos rodar nosso modelo de análise de variância.
# Modelo
Modelo1 <- aov(prod ~ rep + prog,
contrasts = list(prog = "contr.sum"),
data = Colheita_1)
anova(Modelo1)## Analysis of Variance Table
##
## Response: prod
## Df Sum Sq Mean Sq F value Pr(>F)
## rep 3 58.90 19.633 2.1836 0.113071
## prog 9 410.32 45.591 5.0708 0.000475 ***
## Residuals 27 242.75 8.991
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Essa análise variância exige alguns pressupostos, podemos verificar eles nos nossos dados usando:
####################################################
###verificar Pressupostos da análise de variância###
####################################################
names(Modelo1)## [1] "coefficients" "residuals" "effects" "rank"
## [5] "fitted.values" "assign" "qr" "df.residual"
## [9] "contrasts" "xlevels" "call" "terms"
## [13] "model"Modelo1_residuals <- Modelo1$residuals #armazenando os erros ou resíduos
# teste de Normalidade DOS ERROS##
#---------------------------------#
shapiro.test (Modelo1_residuals) # Hipótese de Nulidade##
## Shapiro-Wilk normality test
##
## data: Modelo1_residuals
## W = 0.96494, p-value = 0.2461# a hipótese de que os erros são normais, nesse caso, como o p-value = 0.24
# ou seja, é maior que >0.05 ou 0.01 (alfa adotado), não se rejeita a hipotese de normalidade Vamos trabalhar um poquinho com as informações da ANOVA? Primeiro, guardaremos o valor do quadrado médio:
QME <- anova(Modelo1)["Residuals", "Mean Sq"]
QME## [1] 8.990907E a média da primeira colheita da nossa variável resposta (produção):
med <- mean(Colheita_1$prod, na.rm = TRUE)
med## [1] 8.54525Com eles podemos calcular o coeficiente de variação (CV):
CVe <- (sqrt(QME)/med)*100
CVe## [1] 35.08948Calcule o CVe e QME para a colheita 2
Crie uma função calcular o CVe
Possibilidade de respostas:
CV_E <- function(anova, med){
QME <- anova(anova)["Residuals", "Mean Sq"]
CVe <- (sqrt(QME)/med)*100
return(CVe)
}
##
CV_E(anova = Modelo1, med = med)## [1] 35.08948Podemos também calcular a herdabilidade da característica produção:
n_rep <- nlevels(Colheita_1$rep)
VG <- (anova(Modelo1)["prog", "Mean Sq"]- QME)/n_rep
VE <- QME
H_2 <- VG/ (VG + VE)
H_2## [1] 0.5043871Crie uma função para estimar a herdabilidade
Criando mapas com ggplot
Não vamos entrar em detalhe sobre os códigos que usamos aqui porque ele faz uso do tidyverse, um pacote que permite manipular os dados com muita flexibilidade, mas ele iria requerer outro curso focado apenas nele! Esta é uma aplicação para mostrar que conseguimos fazer quase tudo com R e ggplot em mãos. O mais difícil é saber manipular os dados…
### Instalar pacotes ###
# devtools::install_github("rpradosiqueira/brazilmaps")
# Need to install 'sf'
# install.packages("sf")
# and also the dependencies:
# - 'units' (libudunits2-dev - Linux)
# - 'gdal' (libgdal-dev - Linux)
#
#######################
# Carregar pacote
library(tidyverse)## Warning: package 'tidyverse' was built under R version 4.1.1## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --## v tibble 3.1.2 v dplyr 1.0.7
## v tidyr 1.1.3 v stringr 1.4.0
## v readr 2.0.0 v forcats 0.5.1
## v purrr 0.3.4## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()library(sf)## Linking to GEOS 3.9.0, GDAL 3.2.1, PROJ 7.2.1# Caso não tenha mais conjunto dados disponível no seu ambiente de trabalho
dados = read.csv("dados_alunos2021.csv", stringsAsFactors = FALSE, dec = ",")
# ou
load("dados_alunos2021.RData")
colnames(dados) = c("Ocupacao", "Graduacao", "Conhecimentos_Genetica", "Conhecimentos_Estatistica", "Conhecimento_Gen_Est",
"Latitude", "Longitude")
# Coletando os dados da base por estado
estados = brazilmaps::get_brmap("State")
# Adicionando uma coluna com as siglas dos estados, será importane para a próxima etapa
estados$sigla = c("RO", "AC", "AM", "RR", "PA", "AP", "TO", "MA", "PI", "CE", "RN", "PB", "PE", "AL", "SE", "BA", "MG", "ES", "RJ", "SP", "PR", "SC", "MS", "MT", "GO", "DF", "RS")
## Aqui buscaremos se as coordenadas das cidades estão contidas nas coordenadas
## dos estados
inter_sf = st_as_sf(dados, coords = c('Longitude', 'Latitude'), crs = st_crs(estados))## old-style crs object detected; please recreate object with a recent sf::st_crs()## Mapear coordenadas das cidades nos estados
map_data = st_join(inter_sf, estados %>%
select(sigla, geometry) %>%
mutate(geom_state=geometry)) %>%
as.data.frame() %>%
select(sigla, geom_state)## old-style crs object detected; please recreate object with a recent sf::st_crs()## old-style crs object detected; please recreate object with a recent sf::st_crs()
## old-style crs object detected; please recreate object with a recent sf::st_crs()
## old-style crs object detected; please recreate object with a recent sf::st_crs()
## old-style crs object detected; please recreate object with a recent sf::st_crs()# Verificar quantas pessoas por estado
# Juntar com poligonos dos estados
to_plot = map_data %>%
group_by(sigla) %>%
summarise(n = n()) %>%
right_join(., estados, by=c("sigla"))## old-style crs object detected; please recreate object with a recent sf::st_crs()# Aqui usamos um outro tipo de plot, o geom_sf(), que precisa de dados de
# coordenadas geométricas para construir os mapas. O resto envolve apenas
# funções e argumentos do ggplot!
ggplot(to_plot) +
geom_sf(aes(geometry = geometry, fill = n)) +
theme_minimal() + scale_fill_continuous(low = "#C6FFDD", high = "#f64f59", na.value="white") +
labs(fill = "Número \nde pessoas")
Que tal tentarmos fazer com o mapa mundi?
# Pacotes
# install.packages("rnaturalearth")
# install.packages("rnaturalearthdata")
# install.packages("rgeos")
library("rnaturalearth")
library("rnaturalearthdata")
library("rgeos")
# Mapa mundial
world <- ne_countries(scale = "medium", returnclass = "sf")
## Buscaremos se as coordenadas das cidades estão contidas nas coordenadas
## dos paises
inter_sf = st_as_sf(dados, coords = c('Longitude', 'Latitude'), crs = st_crs(world))
map_data = st_join(world, inter_sf) %>%
filter(!is.na(Graduacao)) %>%
select(name, geometry) %>%
as.data.frame()
# Verificar quantas pessoas por país
# Juntar com poligonos dos países
to_plot = map_data %>%
group_by(name) %>%
summarise(n = n()) %>%
right_join(., world, by=c("name"))
# Aqui usamos um outro tipo de plot, o geom_sf(), que precisa de dados de
# coordenadas geométricas para construir os mapas. O resto envolve apenas
# funções e argumentos do ggplot!
ggplot(to_plot) +
geom_sf(aes(geometry = geometry, fill = n)) +
theme_minimal() + scale_fill_continuous(low = "#C6FFDD", high = "#f64f59", na.value="white") +
labs(fill = "Número \nde pessoas")Família de funções apply
A família de funções apply também podem funcionar como um estrutura de repetição. Sua sintaxe é mais enxuta quando comparada com for ou while e pode facilitar a elaboração do código.
Aqui vamos exemplificar o uso de algumas dessas funções.
apply
A função apply é a base de todas as outras funções da família, portanto a compreensão do funcionamento desta é essencial para entender as demais. Se buscar no help da função, ele indicará que os argumentos da função consistem em: apply(X, MARGIN, FUN, …). Sendo X o conjunto de dados em formato de array (incluindo matrix, que consiste num array de dimensão 2), MARGIN será 1 se a ação deverá ser aplicada à linhas, 2 se for aplicada a colunas e c(1,2) se for aplicada a ambas; FUN é a função que indica ação.
Num simples exemplo temos a matrix:
ex_mat <- matrix(seq(0,21,3), nrow = 2)Se quisermos somar os elementos das colunas usamos:
apply(ex_mat, 2, sum)## [1] 3 15 27 39Se quisermos somar os elementos das linhas:
apply(ex_mat, 1, sum)## [1] 36 48Se fossemos utilizar o for para realizar essa tarefa:
# Soma das colunas
for(i in 1:dim(ex_mat)[2]){
print(sum(ex_mat[,i]))
}## [1] 3
## [1] 15
## [1] 27
## [1] 39# Soma das linhas
for(i in 1:dim(ex_mat)[1]){
print(sum(ex_mat[i,]))
}## [1] 36
## [1] 48lapply
Se diferencia do apply por poder receber outros tipos de objetos (mais utilizado com listas) e devolver o resultado em uma lista.
ex_list <- list(A=matrix(seq(0,21,3), nrow = 2),
B=matrix(seq(0,14,2), nrow = 2),
C= matrix(seq(0,39,5), nrow = 2))
str(ex_list)## List of 3
## $ A: num [1:2, 1:4] 0 3 6 9 12 15 18 21
## $ B: num [1:2, 1:4] 0 2 4 6 8 10 12 14
## $ C: num [1:2, 1:4] 0 5 10 15 20 25 30 35Para selecionar a segunda coluna de todas as matrizes
lapply(ex_list, "[", 2)## $A
## [1] 3
##
## $B
## [1] 2
##
## $C
## [1] 5sapply
A função sapply funciona como o lapply a diferença é que ele retorna apenas um valor por componente da lista e os deposita em um vetor de resposta. Como no exemplo:
sapply(ex_list, "[",1,3)## A B C
## 12 8 20tapply
Esta função é um pouco diferente das demais, ela exige que exista alguma variável categórica (fator) para aplicar ação separadamente conforme suas categorias (levels). Por isso, normalmente é aplicada a data.frames.
Vamos utilizar nosso conjunto de dados:
str(dados)## 'data.frame': 116 obs. of 7 variables:
## $ Ocupacao : chr "Doutorado" "Doutorado" "Mestrado" "Outro" ...
## $ Graduacao : chr "Agronomia" "Agronomia" "Biotecnologia" "Licenciatura em Ciências Biológicas" ...
## $ Conhecimentos_Genetica : chr "Avançado" "Avançado" "Intermediário" "Intermediário" ...
## $ Conhecimentos_Estatistica: chr "Intermediário" "Intermediário" "Básico" "Intermediário" ...
## $ Conhecimento_Gen_Est : chr "Intermediário" "Intermediário" "0" "Intermediário" ...
## $ Latitude : num -25.5 -22.7 -15.9 52.1 -22.8 ...
## $ Longitude : num -54.6 -47.6 -48.1 -106.7 -47.7 ...dados$Graduacao <- as.factor(dados$Graduacao)
dados[dados == "Avançado"] <- 3
dados[dados == "Intermediário"] <- 2
dados[dados == "Básico"] <- 1
dados[dados == "Sabe o teste F? Fui eu que desenvolvi, mas não publiquei porque sou modesto"] <- 99
dados$Conhecimentos_Genetica <- as.numeric(dados$Conhecimentos_Genetica)
tapply(dados$Conhecimentos_Genetica, dados$Graduacao, mean)## Agronomia
## 2.318182
## AGRONOMIA
## 2.000000
## Bacharel em Agronomia
## 2.000000
## Biologia
## 2.571429
## Biología
## 2.000000
## Biotecnologia
## 2.600000
## Ciências bilógicas
## 3.000000
## Ciencias Biologicas
## 2.000000
## Ciências Biológicas
## 2.333333
## Eng. Agronômica
## 2.000000
## Eng. Florestal
## 2.000000
## Engeheria em Biotecnologia Vegetal (Chile)
## 3.000000
## Engenharia Agrícola
## 3.000000
## Engenharia Agronômica
## 1.800000
## Engenharia Agronômica / Licenciatura em Ciências Agrárias
## 2.000000
## Engenharia Biotecnológica
## 2.000000
## Engenharia Florestal
## 2.000000
## Engenheiro Agrônomo
## 2.000000
## Engenheiro Florestal
## 1.000000
## FAEM - UFPEL
## 2.000000
## Genética e Melhoramento
## 2.500000
## Genética e Melhoramento de Plantas
## 2.500000
## Licenciatura em Ciências Biológicas
## 2.000000
## Melhoramento genético
## 2.000000
## Zootecnia
## 2.000000Saiba mais sobre essa família de funções no link
Observe que nas funções apply podemos trocar as funções prontas do r base por funções personalizadas.
Se estiver adiantada/o em relação aos colegas, você já pode fazer os exercícios da Sessão extra, se não, faça-os no conforto do seu lar e nos envie dúvidas pelo fórum.
Gerando sorteio para delineamento experimental
Que tal pensarmos em um delineamento experimental e sortearmos as unidades experimentais? Com o pacote agricolae elaboraremos o sorteio de um experimento de blocos completos casualizados.
############################
#SORTEIO DE EXPERIMENTOS####
############################
#install.packages("agricolae")
library(agricolae)
trt <- c("0","1","2","5","10","20","50","100","Dina")
rcbd <- design.rcbd(trt, 6, serie = 1, seed = 1, "default") # seed = 1
rcbd # Planilha de campo## $parameters
## $parameters$design
## [1] "rcbd"
##
## $parameters$trt
## [1] "0" "1" "2" "5" "10" "20" "50" "100" "Dina"
##
## $parameters$r
## [1] 6
##
## $parameters$serie
## [1] 1
##
## $parameters$seed
## [1] 1
##
## $parameters$kinds
## [1] "default"
##
## $parameters[[7]]
## [1] TRUE
##
##
## $sketch
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
## [1,] "Dina" "5" "50" "0" "1" "20" "2" "100" "10"
## [2,] "1" "2" "100" "0" "10" "20" "Dina" "50" "5"
## [3,] "50" "0" "Dina" "10" "20" "100" "5" "1" "2"
## [4,] "10" "1" "100" "20" "0" "5" "2" "Dina" "50"
## [5,] "2" "20" "1" "50" "5" "10" "100" "Dina" "0"
## [6,] "50" "20" "0" "5" "Dina" "2" "1" "100" "10"
##
## $book
## plots block trt
## 1 11 1 Dina
## 2 12 1 5
## 3 13 1 50
## 4 14 1 0
## 5 15 1 1
## 6 16 1 20
## 7 17 1 2
## 8 18 1 100
## 9 19 1 10
## 10 21 2 1
## 11 22 2 2
## 12 23 2 100
## 13 24 2 0
## 14 25 2 10
## 15 26 2 20
## 16 27 2 Dina
## 17 28 2 50
## 18 29 2 5
## 19 31 3 50
## 20 32 3 0
## 21 33 3 Dina
## 22 34 3 10
## 23 35 3 20
## 24 36 3 100
## 25 37 3 5
## 26 38 3 1
## 27 39 3 2
## 28 41 4 10
## 29 42 4 1
## 30 43 4 100
## 31 44 4 20
## 32 45 4 0
## 33 46 4 5
## 34 47 4 2
## 35 48 4 Dina
## 36 49 4 50
## 37 51 5 2
## 38 52 5 20
## 39 53 5 1
## 40 54 5 50
## 41 55 5 5
## 42 56 5 10
## 43 57 5 100
## 44 58 5 Dina
## 45 59 5 0
## 46 61 6 50
## 47 62 6 20
## 48 63 6 0
## 49 64 6 5
## 50 65 6 Dina
## 51 66 6 2
## 52 67 6 1
## 53 68 6 100
## 54 69 6 10Podemos exportar e salvar nosso sorteio com:
write.table(rcbd,"SORTEIO.txt", row.names=FALSE, sep="\t")
file.show("SORTEIO.txt")
write.csv(rcbd,"SORTEIO.csv",row.names=F)PCA - Análise de componentes principais
O PCA é uma análise exploratória muito frequente quando utilizada por usuários que possuem muitos dados. Nela, torna-se possível visualizar dados que deveriam ter diversas dimensões em apenas duas componentes mais informativas. Para isto, usamos a função autoplot().
# para instalar, use install.packages("ggfortify")
library(ggfortify)
autoplot(prcomp(~ conhecimentoR + altura + peso + idade, data = dados),
data = dados, shape = 'area', colour = 'Graduacao', loadings = TRUE, loadings.colour = 'red',
loadings.label = TRUE, loadings.label.colour = 'black',
loadings.label.size = 4) +
labs(shape = "Área de \ninteresse", color = "Formação", title = "PCA")
# Os parâmetros relacionados com "LOADINGS" são na verdade a direção onde a váriável mais cresce em duas dimensões
# Ela serve apenas como uma aproximação, já que os dados reais precisariam de mais dimensões, que nosso olho não pode verAgora que você fez um plot de PCA, você pode perceber que ele é apenas uma variação do geom_point(), isto permite que você adicione todas as camadas de ggplot que você quiser sobre este gráfico. Já fizemos isto com as legendas!
Sugestões, críticas e elogios
Caso tenha sugestões para aprimoramento desse material, enviar e-mail para biometriamarcadores@gmail.com.
Acesse também outros materiais em português produzidos pela mesma equipe aqui.
Este material foi produzido por alunos do programa de pós-graduação em Genética e Melhoramento de Plantas. Cristiane Taniguti, Fernando Correr e e Wellingson Araújo ministraram o Treinamento.
Também recomendamos materiais em inglês aqui.