Ajustando a distribuição t no parâmetro R: scaling

17

Como ajusto os parâmetros de uma distribuição t, ou seja, os parâmetros correspondentes à 'média' e 'desvio padrão' de uma distribuição normal. Eu suponho que eles são chamados de 'média' e 'escala / graus de liberdade' para uma distribuição t?

O código a seguir geralmente resulta em erros de "falha na otimização".

library(MASS)
fitdistr(x, "t")

Preciso escalar x primeiro ou converter em probabilidades? Qual a melhor forma de fazer isso?

user12719
fonte
2
Ele falha não porque você precisa dimensionar parâmetros, mas porque o otimizador falha. Veja minha resposta abaixo.
Sergey Bushmanov

Respostas:

16

fitdistrusa técnicas de máxima probabilidade e otimização para encontrar parâmetros de uma determinada distribuição. Às vezes, especialmente para a distribuição t, como @ user12719 notou, a otimização no formato:

fitdistr(x, "t")

falha com um erro.

Nesse caso, você deve ajudar o otimizador, fornecendo ponto de partida e limite inferior para começar a procurar por parâmetros ideais:

fitdistr(x, "t", start = list(m=mean(x),s=sd(x), df=3), lower=c(-1, 0.001,1))

Observe, df=3é o seu melhor palpite sobre o que dfpoderia ser um "ótimo" . Depois de fornecer essas informações adicionais, seu erro desaparecerá.

Alguns trechos para ajudá-lo a entender melhor a mecânica interna de fitdistr:

Para as distribuições Normal, log-Normal, geométrica, exponencial e Poisson, os MLEs de forma fechada (e erros padrão exatos) são usados ​​e startnão devem ser fornecidos.

...

Para as seguintes distribuições nomeadas, os valores iniciais razoáveis ​​serão computados se startomitidos ou apenas parcialmente especificados: "cauchy", "gama", "logistic", "binomial negativo" (parametrizado por mu e tamanho), "t" e "weibull " Observe que esses valores iniciais podem não ser bons o suficiente se o ajuste for ruim: em particular, eles não são resistentes a valores discrepantes, a menos que a distribuição ajustada seja de cauda longa.

Sergey Bushmanov
fonte
1
Ambas as respostas (Flom e Bushmanov) são úteis. Estou escolhendo este, porque torna mais explícito que com os valores iniciais corretos e as restrições a otimização 'fitdistr' converge.
user12719
10

νt

νt

set.seed(1234)
n <- 10
x <- rt(n,  df=2.5)

make_loglik  <-  function(x)
    Vectorize( function(nu) sum(dt(x, df=nu,  log=TRUE)) )

loglik  <-  make_loglik(x)
plot(loglik,  from=1,  to=100,  main="loglikelihood function for df     parameter", xlab="degrees of freedom")
abline(v=2.5,  col="red2")

insira a descrição da imagem aqui

n

Vamos tentar algumas simulações:

t_nu_mle  <-  function(x) {
    loglik  <-  make_loglik(x)
    res  <-  optimize(loglik, interval=c(0.01, 200), maximum=TRUE)$maximum
    res   
}

nus  <-  replicate(1000, {x <- rt(10, df=2.5)
    t_nu_mle(x) }, simplify=TRUE)

> mean(nus)
[1] 45.20767
> sd(nus)
[1] 78.77813

Mostrar que a estimativa é muito instável (olhando para o histograma, uma parte considerável dos valores estimados está no limite superior dado para otimizar 200).

Repetindo com um tamanho de amostra maior:

nus  <-  replicate(1000, {x <- rt(50, df=2.5)
    t_nu_mle(x) }, simplify=TRUE)
> mean(nus)
[1] 4.342724
> sd(nus)
[1] 14.40137

o que é muito melhor, mas a média ainda está muito acima do valor real de 2,5.

Lembre-se de que esta é uma versão simplificada do problema real, onde os parâmetros de localização e escala também precisam ser estimados.

tν

kjetil b halvorsen
fonte
5
Sua conclusão de que os problemas de estimativa de df podem realmente funcionar contra o motivo de escolher uma distribuição t em primeiro lugar (ou seja, robustez) é instigante.
User12719
1
(+1) "Não vinculado acima" não é uma resposta errada e pode ser útil para alguns propósitos quando associado a uma estimativa de intervalo. O importante é não usar cegamente as informações de Fisher observadas para formar intervalos de confiança de Wald.
Scortchi - Restabelece Monica
8

Na ajuda para fitdistr, este exemplo:

fitdistr(x2, "t", df = 9)

indicando que você só precisa de um valor para df. Mas isso pressupõe padronização.

Para mais controle, eles também mostram

mydt <- function(x, m, s, df) dt((x-m)/s, df)/s
fitdistr(x2, mydt, list(m = 0, s = 1), df = 9, lower = c(-Inf, 0))

onde os parâmetros seriam m = média, s = desvio padrão, df = graus de liberdade

Peter Flom - Restabelece Monica
fonte
1
Acho que estou confuso sobre os parâmetros de uma distribuição t. Possui 2 (média, df) ou 3 (média, desvio padrão, df)? Eu queria saber se alguém poderia ajustar o parâmetro 'df'.
user12719
1
@ user12719 A distribuição t de Student tem três parâmetros: localização, escala e graus de liberdade. Eles não são referidos como média, desvio padrão e df porque a média e a variação dessa distribuição dependem dos três parâmetros. Além disso, eles não existem em alguns casos. Peter Flom está corrigindo o df, mas isso também pode ser considerado um parâmetro desconhecido.
1
@PeterFlom No caso da distribuição Cauchy , é explícito que m e s são a localização e a escala. Concordo que a notação me sugere que eles representam a média e o desvio padrão, respectivamente. Mas isso pode ser apenas uma simplificação \mue \sigmatambém. +1 a muito tempo atrás, a propósito.
1
@PeterFlom Esta citação do arquivo de ajuda do R implica que df é sempre 9 para distribuição dos alunos? Você não acha que o df também deve ser estimado? Na verdade, a ausência de dfé a causa do erro, e a resposta certa deve fornecer alguma receita para encontrá-lo.
Sergey Bushmanov
1
@ PeterFlom BTW, se você ler o arquivo de ajuda algumas linhas acima da sua citação, descobrirá por que df=9é bom no exemplo deles e irrelevante aqui.
Sergey Bushmanov