<WRAP tabs>
  * [[02_tutoriais:tutorial2:start|Tutorial]]
  * [[01_curso_atual:exercicios2| Exercícios]]
  * [[03_apostila:03-funcoes| Apostila]] 
</WRAP>


====== 2. Tutoriais de Funções Matemáticas no R ======

{{:02_tutoriais:tutorial2:teachMath.gif?400  |}}


No nosso [[02_tutoriais:tutorial1:start|Tutorial 1a]] de introdução à linguagem, vimos alguns tipos básicos de informações que podem ser armazenados em objetos do R e a **estrutura básica de armazenamento que são os vetores**. Neste tutorial vamos apresentar algumas ferramentas de operações algébricas e estatísticas no R.  

<WRAP center round box 60%>

{{youtube>rLjtCxxWR9U}}

</WRAP>

===== Criando Vetores =====

As características básica do vetor  no R são: 

  * uma dimensão de dados 
  * tamanho definido pelo número de elementos
  * armazena dado de um só tipo
  * classificado pelo tipo de dado armazenado 

Os vetores podem ser construídos de várias formas no R, as principais são:
<WRAP center round box 70%>

^função ^ nome ^descrição ^
| <code rsplus> c() </code> | combinar | concatena elementos|
| <code rsplus> seq() </code> | sequência | sequência de valores |
| <code rsplus> : </code> | //colon// | sequência de inteiros |
| <code rsplus> rep() </code> | replicar | repete elementos |


</WRAP>


Execute o código abaixo para entender como criar uma sequência de 1 a 10 em intervalos inteiros de diferentes maneiras:

<code rsplus>
a <- c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10)
b <- seq(from = 1, to = 10, by = 1)
c <- 1:10
length(a)
class(a)
length(c)
class(c)
a == b
a == c
identical(a, b)
identical(a, c)

</code>


Acima, apesar dos objetos conterem a mesma informação de dados, o objeto ''c''  pertence a classe
''integer'', representada por números inteiros, enquanto os outros objetos são da classe ''numeric''. Os objetos armazenam a mesma coisa, mas não são idênticos já que o atributo de classe é diferente. Essa atribuição a classe é feita pelas funções ao criar os objetos.   

===== Operações com Vetores =====

As operações algébricas com vetores apresentam duas características importantes: **equivalência de posição** e **ciclagem do vetor de menor tamanho**. A **equivalência** define que os valores serão operados respeitando suas posições equivalentes entre dois ou mais vetores de mesmo tamanho. No caso de vetores de tamanhos diferentes, há a ciclagem dos elementos do vetor de menor tamanho. Nesse último caso, se o vetor maior é múltiplo do menor, não há nenhuma mensagem na operação. Caso não possa ciclar todos os elementos o mesmo número de vezes, o R opera parte do vetor menor no último ciclo e avisa com uma mensagem de alerta: 'longer object length is not a multiple of shorter object length'. 

Estude os resultados dos comandos nos três exemplos abaixo para entender as operações com vetores.

<code rsplus>
a <- 1:10
b <- -(1:10)
a+b
sum(a, b)
</code>


Aqui temos a característica da **equivalência**. Os vetores de mesmos tamanhos, quando operados algebricamente, operam os valores das posições equivalentes. 


<code rsplus>
a <- seq(from = 0, to = 1, length = 9)
b <- seq(from = 0, to = 0.5, length = 9)
a
b
a/b
b/a
1+b/a
(1+b)/a
</code>

Nesse último bloco de código temos, além da equivalência entre as posições de vetores de mesmo tamanho, algumas operações algébricas que retornam palavras reservadas no R: ''NaN'', ''Inf'' e ''-Inf'':
O ''NaN'' significa 'not a number', ou seja um valor numérico indefinido ou irrepresentável. Operações como ''0/0'' ou ''√-1'' retornam ''NaN''. Dividir um número  real por ''0'' retorna ''Inf'' ou ''-Inf''((apesar de divisão por 0 não ser definida matematicamente, o R considera que o valor a ser dividido é muito próximo de 0)).  


<code rsplus>
a <- rep( c(1, 2), each = 3 )
length(a)
b <- rep( c(3, 4), 6)
length(b)
a+b
</code>


<code rsplus>
a <- rep( c(1, 2, 3), each = 3 )
length(a)
b <- rep( c(3, 4), length = 7)
length(b)
a+b
</code>

Nestes blocos temos o princípio da ciclagem. Ao realizar operações com vetores de tamanhos diferentes, o vetor menor é ciclado, ou seja, reutilizado para a operação. Note que se o tamanho dos vetores for múltiplo, todas as posições do vetor serão utilizadas um mesmo número de vezes e o R não retorna nenhuma mensagem de aviso((lembre-se de sempre inspecionar os objetos que você criar!)). Se os vetores não forem múltiplos a operação é realizada porém com uma mensagem de aviso.


===== Operações sintéticas =====

As operações de vetores acima retornam vetores do tamanho do vetor operado, ou do maior vetor.  Vamos ver algumas operações que retornam vetores de tamanhos distintos ou apenas um atributo do vetor. Para isso, vamos utilizar a seguinte informação de uma pessoa em dieta que anotou seu peso mensalmente durante um ano e obteve os valores:
<code rsplus>
pesos <- c(78.4, 79.8, 76.0, 75.3, 77.4, 78.6, 77.9, 78.8, 79.2, 75.2, 75.0, 79.4)
</code>

A diferença de peso entre um mês e o consecutivo é obtida pela função ''diff'':
<code rsplus>

pesos.dif <- diff(pesos)

</code>

Que tem a particularidade de retornar como resultado um vetor de tamanho igual ao comprimento do vetor de entrada, menos uma posição.

<code rsplus>
length(pesos)
length(pesos.dif)
</code>

Vamos calcular o valor mínimo, máximo e médio do peso:

<code rsplus>
max(pesos)
min(pesos)
mean(pesos)
</code>

O valor mínimo e máximo também é obtido com:
<code rsplus>
range(pesos)
</code>

E mesmo que não tivéssemos a função ''mean'' poderíamos calcular a média com:
<code rsplus>
sum(pesos)/length(pesos)
</code>


<WRAP center round tip 90%>

===== Vocabulário de Referência =====

[[http://cran.r-project.org/doc/contrib/Short-refcard.pdf|{{:02_tutoriais:tutorial2:cardR.png?800|  }}]]

Um dos problemas que enfrentamos ao aprender uma linguagem é a falta de vocabulário que nos impede de fazer uma comunicação eficiente. Isso acontece no R também. Além disso, é mais fácil encontrar um documento que explica como fazer uma análise complexa do que encontrar o nome de uma função que faz uma tarefa simples! É muito bom poder contar com um dicionário da linguagem nessas horas. No caso do R, e outras linguagens computacionais, é comum termos cartões de referência da linguagem. São, geralmente, poucas páginas contendo um vocabulário básico para se comunicar. Na imagem temos a página com as referências às funções matemáticas. Pode baixar o arquivo clicando na imagem do card ou no nosso [[01_curso_atual:material:start#rcard|material de apoio]].
**__É possível fazer quase tudo que precisamos no R com essas 4 páginas de vocabulário básico do card!!__** Imprima frente e verso e mantenha junto ao computador!

</WRAP>
 

===== Operando NaN =====

Calcule a média do logaritmo na base 10 das diferenças de peso, obtidas no tópico anterior ([[02_tutoriais:tutorial2:start#operações_sintéticas|"Operações sintéticas"]]): 

<code rsplus>
mean(log(pesos.dif, base = 10))
</code>

Este comando retorna um aviso e um resultado não numérico, ''NaN'', pois não existem logaritmos de números negativos((pois retorna números complexos)). ''NaN'' significa //Not a Number//, e serve para alertar que o usuário tentou realizar uma operação matemática não definida numericamente ou valores não pertencentes aos números reais.

<code rsplus>
pesos.dif
log(pesos.dif, base = 10)
</code>

===== Operando Valores Faltantes =====
 
É muito comum termos posições em vetores que não tem informação disponível. Por exemplo, no caso acima a pessoa pode ter esquecido de anotar o peso em um dos meses. Nesses casos, para registrar que o valor não está disponível podemos usar a palavra reservada ''NA''.
Basta um valor faltante (''NA'', que significa //Not Available//) ou não numérico (''NaN'') em um vetor para que operações numéricas retornem ''NA'' e ''NaN'', respectivamente. Verifique:

<code rsplus>
faltaUm <- pesos
faltaUm[5] <- NA
faltaUm 
mean(faltaUm)
</code>

<code rsplus>
sqrt(-1:9)
mean(sqrt(-1:9))
</code>

Muitas funções têm um argumento lógico ''na.rm'' que, se declarado verdadeiro, desconsidera os valores ''NA'' e ''NaN''. Verifique:

<code rsplus>
mean(faltaUm, na.rm = TRUE)
var(sqrt(-1:9), na.rm = TRUE)
</code>


===== Índice de Shannon =====

O índice de diversidade de Shannon é dado pela fórmula:

$$H´=-\sum \ p_i \ ln p_i$$


onde $p_i$ é a proporção de indivíduos da espécie $i$ em relação ao total de indivíduos. A seguir construímos um objeto chamado ''abund'', que representa as abundâncias das espécies em uma comunidade. Em seguida, calculamos o valor do índice, usando objetos intermediários, sempre verificando cada passo intermediário para garantir que os objetos formados contém o valor que queremos:

<code rsplus>
abund <- c(13, 23, 29, 29, 30, 48, 72, 76, 83, 84, 86, 97, 102, 229, 343)
tot <- sum(abund)
tot
pi <- abund/tot
pi
log.pi <- log(pi)
log.pi
pi.log.pi <- pi*log.pi
pi.log.pi
H <- - sum(pi.log.pi)
H
</code>

O mesmo cálculo pode ser feito em uma única linha com:
<code rsplus>
-sum(abund/sum(abund) * log(abund/sum(abund)))
</code>

Compare os resultados.


<WRAP center round box 90%>

===== Funções Aninhadas =====
{{:02_tutoriais:tutorial2:eagle.png?250  |}}

No código acima,  podemos observar uma característica muito básica e fundamental da sintaxe da linguagem R: o aninhamento de funções. O interpretador do R foi desenhado para ler a linha de comando começando pelas funções mais internas e expandindo o resultado para as mais externas, de maneira análoga com o que fizemos ao calcular o valor com os objetos intermediários. No caso do código aninhado, os valores intermediários só existem durante o processamento do código.  


</WRAP>
 

===== Onde foi parar o pi? =====

No [[:02_tutoriais:tutorial2:start#indice_de_shannon|tópico anterior]] criamos um objeto ''pi'', mas este é nome que o R usa para armazenar o número π(([[http://en.wikipedia.org/wiki/Pi|razão entre o perímetro e o diâmetro do círculo]]))!!! 8-O

Mas não se preocupe. O objeto ''pi'' original pertence ao pacote ''base'' e continua lá. 

Como o R busca qualquer objeto chamado na ordem estabelecida no caminho de busca (''search path'') e como o sua área de trabalho (''workspace'') está na frente do pacote ''base'', ele irá usar o novo objeto ''pi'' se você digitar apenas o nome do objeto, pois este será o primeiro objeto com este nome que ele encontrará.
Verifique o caminho de busca (''search path'') com:
<code rsplus>
search()
</code>

Sua área de trabalho (''workspace'') tem o nome padrão  ''.GlobalEnv''  e tem a precedência no caminho de busca de todos os outros compartimentos de memória (''environments''). Agora  o objeto ''pi'' de seu  ''workspace'' tem precedência ao objeto ''pi'' do pacote base.
<code rsplus>
pi
</code>

Podemos indicar para o R que queremos um objeto em um compartimento de memória específico e, nesse caso, ele irá buscar o objeto diretamente. Podemos então acessar o número π explicitando o ambiente (pacote) onde ele está da seguinte forma, ''compartimento::objeto'', como a seguir:

<code rsplus>
base::pi
</code>

Nesse caso, se o objeto não existir nesse compartimento, mas estiver em outro, o R não irá considerá-lo e retornará a mensagem de erro que o objeto não foi encontrado. Veja o código a seguir:

 
<code rsplus>
a <- 0
base::a
</code>


===== Distribuições Probabilísticas =====
<WRAP center round box 60%>

{{youtube>FC-tV1Nlwb8}}

</WRAP>

O R contém muitas distribuições probabilísticas no pacote ''stats'', carregado automaticamente quando iniciamos uma sessão no R. 

Vamos entender como essas distribuições são utilizadas no R com um exemplo de dados.

==== Eita avião apertado ====

{{:02_tutoriais:tutorial2:airplaneseats.jpg?400  |}}Um artigo reportou a largura do quadril de mulheres de cerca de vinte anos como sendo ((Victoria J. Simpson, Gayle Brewer, Colin A. Hendrie. Evidence to Suggest that Women’s Sexual Behavior is Influenced by Hip Width Rather than Waist-to-Hip Ratio. Archives of Sexual Behavior, 2014; DOI: 10.1007/s10508-014-0289-z))  35 ± 3.4 cm (média e desvio padrão). Por outro lado, os assentos de aviões por padrão tem cerca de  40.6 cm ((16 polegadas)) de largura. Com essas informações e algumas premissas podemos inferir muitas coisas. As premissas mais importantes são:

  * a média e o desvio são uma boa estimativa dos parâmetros da distribuição da largura de quadril da população de mulheres de 20 e poucos anos;
  * a distribuição de largura do quadril é bem representada por uma distribuição probabilística normal

A partir desses dois parâmetros (média e desvio) da distribuição normal podemos, por exemplo, simular os dados de uma amostra de 200 tamanhos de quadris no R. A função ''rnorm'' faz a simulação de amostras aleatórias de uma distribuição normal e tem como argumentos ''n'' que é o tamanho amostral, ''mean'' e ''sd'', que são os parâmetros que definem a distribuição normal.

<code rsplus>

amostra <- rnorm(n = 200, mean = 35, sd = 3.4)
hist(amostra)


</code> 

Quantas mulheres na nossa amostra não caberiam nas cadeiras do avião?

<code rsplus>

sum(amostra > 40.6)

</code> 

Conseguiu entender o que foi feito acima? Uma forma de compreender códigos mais complexos é dividir em partes e estudar cada elemento independentemente:


<code rsplus>
vf_41 <-  amostra > 41
vf_41
sum(vf_41) # lembre-se que TRUE pode ser considerado como 1 e FALSE como 0.
</code> 

Como entendemos todos os comandos acima, podemos refazer a simulação de amostra e contar as mulheres que não caberiam nos assentos de avião em uma única linha. Agora vamos repetir a nossa amostra 10x:


<code rsplus>

sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)
sum(rnorm(n = 200, mean = 35, sd = 3.4) > 40.6)

</code> 


A cada amostra nosso valor varia um pouco. Isso ocorre porque a simulação de amostra é uma realização aleatória de 200 valores da distribuição normal teórica. Podemos obter o valor teórico dessa proporção de mulheres em nossa população estatística teórica. Para isso usamos a função ''pnorm'' que tem os parâmetros ''q'' de quantil e os parâmetros da normal:

<code rsplus>

pnorm(q = 40.6, mean = 35, sd = 3.4)

</code>

Esse valor é a proporção de mulheres que caberiam nos assentos com 40.6cm de largura. Para o cálculo das que não caberiam podemos fazer de duas formas:

<code rsplus>

1 - pnorm(q = 40.6, mean = 35, sd = 3.4)
pnorm(q = 40.6, mean = 35, sd = 3.4, lower.tail = FALSE)

</code>

Por volta de 5% das mulheres não caberiam nessas poltronas, o que não parece lá muito simpático por parte das companhias aéreas!  
Vamos nos perguntar outra coisa, por exemplo, qual o tamanho de assento no qual caberiam 99% das mulheres de cerca de vinte anos? Para extrair esse valor diretamente da distribuição normal teórica, usamos a função ''qnorm''.

<code rsplus>

qnorm(p = 0.99, mean = 35, sd = 3.4)

</code>

Isso significa que apenas 2.3 cm de aumento na largura das poltronas das aeronaves é capaz de reduzir em 5x a proporção de mulheres que não caberiam. Parece uma boa política para a companhia aérea... 

==== Amplitude de uma Amostra Normal ====

Vamos simular uma amostra de 10 valores, tomados de uma distribuição normal com  média 10 e desvio padrão 2.5:

<code rsplus>
normal.1 <- rnorm(10, mean = 10, sd = 2.5)
</code>

Calcule mínimo e máximo e a amplitude destes valores:
<code rsplus>
range(normal.1)
diff(range(normal.1))
</code>

O que acontece à medida que aumentamos o tamanho da amostra? Verifique para simulações de 100, 1000 e 10000 valores:
<code rsplus>
diff(range(rnorm(100, mean = 10, sd = 2.5)))
diff(range(rnorm(1000, mean = 10, sd = 2.5)))
diff(range( rnorm(10000, mean = 10, sd = 2.5)))
</code> 


==== As Funções que operam distribuições de probabilidades ====


<WRAP center round box 80%>
**If data fails the Teacher's t test, you can just force it to take the test again until it passes.**

[[https://www.explainxkcd.com/wiki/index.php/1347:_t_Distribution|{{  :02_tutoriais:tutorial2:t_distribution.png?600  |}}]]
</WRAP>


O que foi apresentado para **Distribuição Normal**  pode ser generalizado para outras distribuições de probabilidade.

Há quatro principais funções para se extrair informações básicas de distribuições probabilísticas:

  * **d**//distrib// - retorna a //densidade probabilística// para um dado valor da variável;
  * **p**//distrib// - retorna a //probabilidade acumulada// para um dado valor da variável;
  * **q**//distrib// - retorna o //quantil// para um dado valor de probabilidade acumulada;
  * **r**//distrib// - retorna //valores// (números aleatórios) gerados a partir da distribuição;

No caso da Distribuição Normal: //distrib// = ''norm''.  Para outras distribuições temos:

<WRAP center round box 95%>

/* 

<box 70% orange centered | >
*/

**DISTRIBUIÇÕES ESTATÍSTICAS NO R**
^ Distribuição ^ Nome no R ^ Parâmetros((os argumentos de cada função incluem estes parâmetros, entre outras coisas)) ^
| beta 	       | beta      | shape1, shape2, ncp |
| binomial     | binom     | size, prob |
| Cauchy       | cauchy    | location, scale |
| qui-quadrado | chisq     | df, ncp |
| exponential  | exp       | rate |
| F            | f         | df1, df2, ncp |
| gamma        | gamma     | shape, scale |
| geométrica   | geom      | prob |
| hypergeométrica | hyper  | m, n, k |
| log-normal   | lnorm     | meanlog, sdlog |
| logística    | logis     | location, scale |
| binomial negativa | nbinom     | size, prob |
| normal       | norm      | mean, sd |
| Poisson      | pois      | lambda |
| t de Student | t         | df, ncp |
| uniforme     | unif      | min, max |
| Weibull      | weibull   | shape, scale |
| Wilcoxon     | wilcox    | m, n  |


</WRAP>


===== Qui-quadrado na unha =====
Vamos fechar com um exemplo hipotético de um estudo de preferência alimentar. 
Nosso ecólogo virtual estimou a proporção de cinco tipos de itens alimentares para uma espécie de ave em uma área. Esses itens estavam disponíveis na proporção de 60%, 28%, 9%, 2,5% e 0,5%. No mesmo local, amostrou-se ao acaso eventos de alimentação desta ave, contando quantos eventos foram de consumo de cada um dos itens. Os resultados das contagens dos eventos de alimentação foram 544, 285, 117, 54, 12, para cada um dos itens respectivamente.  

Vamos criar os objetos com estes valores: 

<code rsplus>
disp <- c(60, 28, 9, 2.5, 0.5)
cons <- c(544, 285, 117, 54, 12)
</code>

A pergunta que esses trabalhos de dieta se fazem é: "a espécie tem preferência por algum item?". Se isso não acontecer, espera-se que 60% do total de itens consumidos sejam do primeiro tipo, e assim por diante. O total de itens consumidos é:
<code rsplus>
totItens <- sum(cons)
</code>
E os valores esperados pela hipótese de falta de preferência serão:
<code rsplus>
espe <- totItens*disp/100
</code>
O que resulta em uma série de desvios entre os valores esperados e os observados: 
<code rsplus>
desv <- cons - espe
desv
</code>
Para testar esta diferença, calculamos o valor do Qui-quadrado, que é a somatória de cada desvio elevado ao quadrado e dividido pelo respectivo valor esperado:
<code rsplus>
dQuad <- desv^2/espe
dQuad
qui2 <- sum(dQuad)
qui2
</code>
Qual a chance de um valor de Qui-quadrado maior ou igual a este ocorrer por acaso (ou seja, mesmo que não haja preferência)? Como são cinco itens alimentares, temos quatro graus de liberdade, e o que queremos saber é: qual a probabilidade de encontrarmos a diferença observada em relação ao esperado ou maiores, em um cenário onde a espécie não tem preferência alimentar (hipótese nula).   Obtemos isto com a função de probabilidade acumulada da distribuição de Qui-quadrado:
<code rsplus>
pchisq(q = qui2, df = 4, lower.tail = FALSE)
</code>

A probabilidade é baixíssima, o que indica que o consumo não foi na proporção da disponibilidade dos itens alimentares. Usamos o argumento ''lower.tail = FALSE'' para que a função retorne a probabilidade da parte que resta da distribuição, após este valor((mais rigorosamente, da área sob a curva a partir deste valor, isto é, a integral da função de densidade probabilística de Qui-quadrado deste valor até +∞)), que está representado em vermelho na figura abaixo:


{{:02_tutoriais:tut2fig1.png?|}}


== Extra ==
Para reproduzir a figura acima digite os comandos abaixo. Para entender, veja a ajuda da função ''curve''.

<code rsplus>
## Faz o grafico da funcao Qui-quadrado com 4 graus de liberdade,
## veja a ajuda da funcao curve
curve(dchisq(x, df = 4), 0, 70, xlab = "Qui-quadrado, 4 g.l.", ylab = "Densidade probabilística")
## Sobrepoe uma linha vermelha a partir 
##do valor calculado do Qui-quadrado
curve(dchisq(x, df = 4), 56.93, 70, add = T, col = "red", lwd = 2)
</code>

<WRAP center round warning 100%>

<WRAP center round box 90%>


==== Código Canalizado ====

[[https://www.monkeyuser.com/2019/code-entropy/|{{:02_tutoriais:tutorial2:119-code-entropy.png?400  |}}]]

As linguagens de programação, assim como as linguagens naturais, estão em constante modificação. Existe na comunidade do R um movimento com um novo //dialeto// que foi chamado de  ''tidyverse''. Esse conjunto de pacotes propõem formas mais compactas de códigos para manipular e grafar dados((Como toda boa polêmica existem muitos apoiadores e muitos críticos dessa nova sintaxe do R. Para saber mais, veja o artigo [[https://github.com/matloff/TidyverseSkeptic|Tidyverse Skeptic]] de Norman Matloff, professor da Universidade da California. )). Este curso é baseado na forma mais original da linguagem contida nos pacotes básicos da distribuição do R e sem a pretensão de ensinar esse novo dialeto. Entretanto, irão encontrar muita documentação que usa essa filosofia. Uma das dificuldades para entender o código nesse //dialeto// é a introdução do conceito de canalização  de procedimentos, onde o resultado de uma função pode ser direcionada a outra através do //pipe// (''%>%''). Essa ferramenta, que está no pacote ''magrittr'', modifica bastante a lógica de funções aninhadas do interpretador do R. Abaixo listamos alguns dos conceitos básicos para entender os códigos em ''tidyverse'':

  * ''x %>% f'' é equivalente a ''f(x)''
  * ''x %>% f(y)'' é equivalente a ''f(x, y)''
  * ''x %>% f %>% g %>% h '' é equivalente a ''h(g(f(x)))''
  * ''x %>% f(y, .)'' é equivalente a ''f(y, x)''
  * ''x %>% f(y, z = .)'' é equivalente a ''f(y, z = x)''

 
</WRAP>

</WRAP>