Funções de Densidade Probabilísticas (PDF em inglês) descrevem valores de densidade probabilística que uma função matemática f(x) atribui a cada valor x. No caso de variáveis discretas, como contagens, as PDFs atribuem uma probabilidade a cada valor:
$$f(x) = P(X=x)$$ 1)
Para variáveis contínuas, a PDF atribui densidades de probabilidade a cada valor. A integral dessas densidades em um intervalo dá a probabilidade da variável ter o valor nesse intervalo:
$$\int_a^b f(x) = P(a\leqx \leq b)$$ 2)
A Função de Probabilidade Cumulativa (CDF em inglês), por sua vez, retorna a probabilidade acumulada até o valor x:
$$F(x) = \int_\Omega^x f(x) = P(X\leqx)$$ 3)
Mas como obter uma PDF a partir de uma função f(x) qualquer? Para tanto é necessário transformar a função para que ela retorne probabilidades ao invés de valores na mesma unidade de x. Para isso, a função deve somar um se for uma PDF discreta, ou integrar um, se for contínua. Isso garante uma das propriedades necessárias de uma PDF: a soma da probabilidade de todos os eventos possíveis é um.
Para PDFs contínuas, temos então que dividir a função pela sua integral definida. Assim, podemos obter uma PDF dividindo a função f(x) por ∫ f(x) dx. Dessa forma, a nova função terá integral igual a 1 (um).
Em ecologia, observações são geralmente números positivos. Alguns fenômenos possuem observações com valor zero enquanto outros só são observáveis com valores maiores que zero. Assim, a integral usada para obter a PDF de uma função é geralmente definida entre zero e infinito.
A função integrate do R encontra a integral definida de uma função por integração numérica. Porém, estamos interessados na integral simbólica da função, ou seja, buscamos qual função cuja derivada seja nossa função inicial. Assim, é possível utilizar o software livre Maxima (http://maxima.sourceforge.net/) que trabalha com álgebra e cálculo simbólico.
Neste tutorial iremos usar a função de Ricker como exemplo:
$$f(x) = a x e^{-b*x}$$
Aqui, x é um valor qualquer, enquanto a e b são os parâmetros da função. Trata-se de uma função cujos valores de f(x) crescem rapidamente, atingem um pico e depois decaem com menor intensidade (Figura 1). Pode ser uma boa função para descrever fenômenos que atingem um pico e depois decaem, como fecundidade em relação a idade de fêmeas, por exemplo.
Figura 1. Exemplo da função de Ricker definida para a= 2 e b= 6.
A Ricker é uma função bastante simples que atende aos objetivos deste tutorial. Funções mais complexas e com mais parâmetros podem não se 'comportar' tão bem, ou podem ter integrais complicadas. Em muitos casos, para encontrar a PDF também será necessário definir o domínio dos parâmetros, pois em muitos casos há soluções distintas para domínios diferentes.
Uma vez obtida a PDF da função, basta integrá-la de zero a x para encontrar a função que descreve a sua CDF. As CDF são úteis caso se pretenda truncar uma PDF, como veremos abaixo, no passo 6.
ricker = function (x,a,b,...) {a*x*exp(-b*x)}
integrate (a*x*exp(-b*x), x, 0, inf) > a/b^2
Esse código tem como argumentos a função (Ricker), a variável de interesse (x), e os valores inferior (zero) e superior (infinito) de integração.
[a*x*exp(-b*x)] / [a/b^2] > b^2*x*exp(-b*x)
Note que o parâmetro a já não existe mais na PDF
dricker = function (x,b,...) {b^2*x*exp(-b*x)} integrate (dricker, lower=0, upper=Inf, b=2)
Varie os valores de b. Se a PDf estiver correta, a integração numérica deve retornar valores muitíssimo próximos de 1.
integrate (b^2*x*exp(-b*x), x, 0, x) > b^2*(1/b^2-((b*x+1)*%e^(-b*x))/b^2) # simplificando... > exp(-b*x)*(exp(b*x)-b*x-1)
pricker = function (x, b,...){exp(-b*x)*(exp(b*x)-b*x-1)} # CDF dricker.t = function (x, b,...){ p.trunc = pricker(trunc,b) dricker (x,b)/(1-p.trunc)}
Aqui usamos a CDF para encontrar a probabilidade acumulada até o ponto de truncagem desejado (i.e.’ trunc’) , para que possamos corrigir a PDF ‘adicionando’ essa probabilidade ao restante da PDF. Um teste para verificar se sua PDF truncada está funcionando é verificar se ele gera valores de probabilidades maiores que a PDF original.
dricker (10, b=0.1) > 0.03678794 dricker.t (10, b=0.1) > 0.04013735 # OK!
An unforeseen error has occured. This is most likely a bug somewhere.
More info has been written to the DokuWiki error log.