Ecological Analyses in Julia

Eduardo Santos

doi:10.5281/zenodo.23126524

4 Introdução ao Julia

Este capítulo apresenta a linguagem para quem já programa em R. O foco está nas diferenças que importam na prática, não em uma varredura exaustiva da sintaxe.

using Statistics
using Random
Random.seed!(42)

TaskLocalRNG()

Como funcionam os exercícios

A maioria das seções termina com um quadro Sua vez. Resolva no VS Code, com a pasta do livro aberta (Seção 3.8), e carregue o verificador uma vez por sessão:

include("exercises/checar.jl")

Depois passe sua resposta — um valor, ou uma função que você escreveu — junto com o número do exercício:

checar("4.3", simpson)

checar diz se a resposta está certa e, se não estiver, o que ela fez e o que era esperado. Cada quadro tem uma dica e uma solução recolhidas: tente primeiro, depois espie.

4.1 Vetores e indexação

alturas = [12.4, 8.1, 15.7, 9.3, 11.0]

length(alturas), first(alturas), last(alturas)

(5, 12.4, 11.0)

Separar valores por vírgulas cria uma tupla: um grupo fixo de valores, mostrado entre parênteses. É assim que Julia devolve vários resultados de uma vez.

Duas coisas funcionam como em R, e uma não:

alturas[1]        # indexação começa em 1, como em R
alturas[2:4]      # intervalos são inclusivos nas duas pontas

3-element Vector{Float64}:
  8.1
 15.7
  9.3

# Julia NÃO recicla vetores silenciosamente.
# try ... catch captura o erro e o mostra, em vez de interromper o livro:
try
    [1, 2, 3] + [1, 2]
catch e
    e
end

DimensionMismatch("a has size (3,), b has size (2,), mismatch at dim 1")

Em R

c(1,2,3) + c(1,2) recicla o vetor mais curto e emite apenas um aviso; Julia lança DimensionMismatch, como o bloco acima demonstra. A diferença é de política: R prioriza seguir adiante, Julia prioriza falhar cedo. Qual serve melhor depende do seu fluxo de trabalho — o que importa é saber qual das duas você está usando.

Sua vez 4.1 — Tire um valor errado

Estas são as áreas (m²) de cinco lagoas. A última é um erro de digitação.

areas = [174.1, 161.6, 233.9, 174.5, 1980.0]

5-element Vector{Float64}:
  174.1
  161.6
  233.9
  174.5
 1980.0

Crie limpo, com todas as lagoas menos a última. Não digite a posição 4: seu código deve continuar funcionando se mais lagoas forem adicionadas.

Dica

Dentro dos colchetes, end significa “a última posição”; logo, end - 1 é a anterior.

Solução

limpo = areas[1:end-1]
checar("4.1", limpo)

✓ Correto!

4.2 Broadcasting: o ponto

Em R, quase toda função é vetorizada por padrão. Em Julia, você declara a vetorização com um ponto. Isso é explícito e vale para qualquer função, inclusive as que você escreve.

log.(alturas)              # aplica log elemento a elemento
alturas .^ 2               # potência elemento a elemento
alturas .> 10              # comparação elemento a elemento

5-element BitVector:
 1
 0
 1
 0
 1

# funciona igual com uma função sua
padroniza(x, μ, σ) = (x - μ) / σ

padroniza.(alturas, mean(alturas), std(alturas))

5-element Vector{Float64}:
  0.3724002529860428
 -1.0833461905048525
  1.4896010119441712
 -0.6770913690655326
 -0.10156370535983013

A regra do ponto

Se você quer “aplicar isso a cada elemento”, ponha um ponto. f.(x), x .+ y, x .== y. Uma cadeia inteira pode ser pontuada de uma vez com @., um macro (Seção 4.8):

@. log(alturas) + 1

5-element Vector{Float64}:
 3.517696472610991
 3.0918640616783932
 3.7536607123542622
 3.2300144001592104
 3.3978952727983707

Sua vez 4.2 — Log das áreas

O gráfico de área em Seção 3.8 usou escala log, porque as áreas das lagoas são muito assimétricas à direita. Calcule o logaritmo na base 10 de cada valor de areas e guarde em log_areas.

Dica

log10(areas) dá erro. Leia a mensagem e releia a regra do ponto logo acima.

Solução

log_areas = log10.(areas)
checar("4.2", log_areas)

✓ Correto!

log10 é definida para um número; o ponto a aplica a cada elemento.

4.3 Funções

# forma curta, para funções de uma linha
riqueza(abundancias) = count(>(0), abundancias)

# forma longa
function shannon(abundancias)
    total = sum(abundancias)
    p = abundancias ./ total
    p = filter(>(0), p)          # log(0) é indefinido
    return -sum(p .* log.(p))
end

comunidade = [45, 12, 0, 8, 23, 0, 3]
riqueza(comunidade), shannon(comunidade)

(5, 1.289250407583511)

>(0) é uma função que pergunta “é maior que 0?”. count(>(0), x) conta quantos elementos de x respondem que sim, e filter(>(0), x) mantém só esses.

Sua vez 4.3 — Índice de Simpson

Escreva simpson(abundancias), que retorna o índice de Simpson $1 - \sum p_i^2$, onde $p_i$ é a proporção de indivíduos da espécie $i$. As primeiras linhas de shannon já calculam as proporções.

Antes de chamar checar, teste simpson([3, 1]) você mesmo: deve dar 0.375.

Dica

p .^ 2 eleva cada proporção ao quadrado e sum soma tudo. Zeros não precisam de cuidado especial aqui: $0^2 = 0$.

Solução

function simpson(abundancias)
    p = abundancias ./ sum(abundancias)
    return 1 - sum(p .^ 2)
end

checar("4.3", simpson)

✓ Correto!

4.4 Despacho múltiplo

Esta é a ideia central da linguagem, e a razão pela qual pacotes Julia compõem sem código de cola. Uma função pode ter vários métodos, e Julia escolhe qual rodar com base nos tipos de todos os argumentos.

# um método para contagens inteiras
diversidade(x::Vector{<:Integer}) = "índice para dados de contagem"

# outro para coberturas proporcionais
diversidade(x::Vector{<:AbstractFloat}) = "índice para dados de cobertura"

diversidade([45, 12, 8]), diversidade([0.45, 0.12, 0.08])

("índice para dados de contagem", "índice para dados de cobertura")

:: declara o tipo que o argumento deve ter, e Vector{<:Integer} quer dizer “um vetor de algum tipo de inteiro” (<: se lê “é um subtipo de”). Um vetor de números decimais não serve para o primeiro método, e Julia escolhe o segundo.

Na prática, isso significa que quando alguém define um novo tipo de dado ecológico, mean, plot e fit passam a funcionar com ele desde que os métodos apropriados sejam definidos — sem que os autores de Statistics ou Makie saibam da existência do novo tipo.

Sua vez 4.4 — Um índice por local

simpson trata uma comunidade. Dados de campo costumam vir como uma matriz locais × espécies, uma linha por local. Adicione um segundo método, simpson(M::Matrix), que retorne um índice por linha. O método para vetores tem que continuar funcionando.

Dica

eachrow(M) percorre as linhas uma a uma. Combine com o ponto, ou com uma compreensão: [f(linha) for linha in eachrow(M)]. Se checar reclamar dos tipos dos argumentos, leia a mensagem com atenção — é uma armadilha comum.

Solução

simpson(M::Matrix) = simpson.(eachrow(M))

locais = [45 12  0 8;
           5  5  5 5;
          30  0  0 0]
simpson(locais)

3-element Vector{Float64}:
 0.47147928994082844
 0.75
 0.0

checar("4.4", simpson)
length(methods(simpson))

✓ Correto!

O método para matrizes reaproveita o método para vetores em cada linha. Uma função, dois métodos, e Julia escolhe pelo tipo do argumento; rode methods(simpson) no REPL para ver os dois.

4.5 Dados faltantes

medidas = [12.4, missing, 15.7, 9.3]

mean(medidas)              # propaga, como em R

missing

mean(skipmissing(medidas))   # equivalente a na.rm = TRUE

12.466666666666669

Em R

NA corresponde a missing. A diferença é que missing é um tipo próprio, e o compilador sabe se uma coluna pode ou não contê-lo — Vector{Float64} nunca contém missing, Vector{Union{Missing,Float64}} pode. Isso elimina uma classe inteira de bugs silenciosos.

Sua vez 4.5 — Quanto está faltando?

A cobertura de dossel não foi medida em todas as lagoas:

dossel = [70.0, 42.6, missing, 49.4, missing, 81.2]

6-element Vector{Union{Missing, Float64}}:
 70.0
 42.6
   missing
 49.4
   missing
 81.2

Escreva resumo_faltantes(x), que retorna uma tupla com dois números: quantos valores estão faltando e a média dos demais. Para dossel, o resultado é (2, 60.8).

Dica

count(f, x) conta os elementos para os quais f é true, e ismissing é uma função desse tipo. A média está no bloco logo acima. Retorne os dois com return (n, m).

Solução

resumo_faltantes(x) = (count(ismissing, x), mean(skipmissing(x)))

resumo_faltantes(dossel)

(2, 60.8)

checar("4.5", resumo_faltantes)

✓ Correto!

4.6 Encadeamento

# o pipe nativo |> encadeia funções de um argumento
comunidade |> shannon |> exp    # número de Hill de ordem 1

3.63006446685595

Para pipelines com múltiplos argumentos, @chain (do Chain.jl, reexportado por TidierData.jl) cumpre o papel de |> do R. Ele é o eixo do Capítulo 5.

4.7 Estruturas de controle

resultado = String[]

for (i, a) in enumerate(comunidade)
    if a == 0
        push!(resultado, "espécie $i ausente")
    elseif a < 10
        push!(resultado, "espécie $i rara")
    else
        push!(resultado, "espécie $i comum")
    end
end

resultado

7-element Vector{String}:
 "espécie 1 comum"
 "espécie 2 comum"
 "espécie 3 ausente"
 "espécie 4 rara"
 "espécie 5 comum"
 "espécie 6 ausente"
 "espécie 7 rara"

String[] cria um vetor de textos vazio, que push! vai enchendo. enumerate entrega a posição i e o valor a de cada elemento. Dentro de um texto, $i insere o valor de i (é a interpolação, como paste0("espécie ", i) em R).

Loops são rápidos

Ao contrário do R, você não precisa vetorizar por performance. Escreva o loop que expressa a ideia. Isso muda como se escrevem modelos nulos e simulações — assunto que volta no Capítulo 11 (veja a Seção 11.6).

Sua vez 4.6 — Shannon, uma espécie por vez

Reescreva shannon como shannon_laco, usando um loop for em vez de broadcasting. O resultado tem que ser igual ao de shannon, inclusive para comunidades com contagens zero, como comunidade.

Dica

Comece com H = 0.0, percorra as contagens uma a uma, transforme cada uma em proporção e subtraia p * log(p) de H. Se aparecer NaN, digite 0 * log(0) no REPL.

Solução

function shannon_laco(abundancias)
    total = sum(abundancias)
    H = 0.0
    for a in abundancias
        a == 0 && continue       # espécies ausentes não contribuem
        p = a / total
        H -= p * log(p)
    end
    return H
end

shannon_laco(comunidade), shannon(comunidade)

(1.289250407583511, 1.289250407583511)

checar("4.6", shannon_laco)

✓ Correto!

a == 0 && continue é a forma curta, em Julia, de if a == 0 continue end.

4.8 Convenções de nomes

Convenção	Significado
`f!(x)`	modifica `x` no lugar (`sort!`, `push!`)
`f.(x)`	aplica `f` elemento a elemento
`@f x`	macro, não função — reescreve código antes de rodar
`Tipo`	tipos em `CamelCase`
`nome_funcao`	funções em minúsculas; sublinhados são bem-vindos quando ajudam a leitura (`resumo_faltantes`)

Desafios

Estes não têm verificador nem uma única resposta certa. São sobre o que o seu código faz — e o que ele não diz.

simpson([-3, 5]) retorna um número. Qual, e deveria? Altere o método de simpson para vetores de modo que contagens negativas gerem um erro informativo com throw(ArgumentError("...")). O método para matrizes herda a correção?
As comunidades [10, 4, 3, 3] e [8, 7, 4, 1] têm a mesma riqueza e o mesmo número de indivíduos. Calcule Shannon e Simpson para as duas. Os dois índices concordam sobre qual comunidade é mais diversa? O que isso implica para um artigo que reporta “diversidade” como um único número?
Carregue data/ponds.csv como em Seção 3.8 e rode resumo_faltantes em ponds.canopy. A contagem bate com o que describe(ponds) mostra?