Ajustement de la distribution t dans R: paramètre d'échelle

17

Comment ajuster les paramètres d'une distribution t, c'est-à-dire les paramètres correspondant à la «moyenne» et à «l'écart-type» d'une distribution normale. Je suppose qu'ils sont appelés «moyenne» et «échelle / degrés de liberté» pour une distribution t?

Le code suivant entraîne souvent des erreurs «échec de l'optimisation».

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

Dois-je d'abord mettre à l'échelle x ou convertir en probabilités? Comment le faire au mieux?

user12719
la source
2
Il échoue non pas parce que vous devez mettre à l'échelle les paramètres, mais parce que l'optimiseur échoue. Voir ma réponse ci-dessous.
Sergey Bushmanov

Réponses:

16

fitdistrutilise des techniques de maximisation de vraisemblance et d'optimisation pour trouver les paramètres d'une distribution donnée. Parfois, en particulier pour la distribution t, comme l'a remarqué @ user12719, l'optimisation sous la forme:

fitdistr(x, "t")

échoue avec une erreur.

Dans ce cas, vous devez donner un coup de main à l'optimiseur en fournissant le point de départ et la borne inférieure pour commencer à rechercher les paramètres optimaux:

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

Notez que df=3c'est votre meilleure estimation de ce que dfpourrait être un "optimal" . Après avoir fourni ces informations supplémentaires, votre erreur disparaîtra.

Quelques extraits pour vous aider à mieux comprendre les mécanismes internes de fitdistr:

Pour les distributions normale, log-normale, géométrique, exponentielle et de Poisson, les MLE de forme fermée (et les erreurs standard exactes) sont utilisées et startne doivent pas être fournies.

...

Pour les distributions nommées suivantes, des valeurs de départ raisonnables seront calculées si elles startsont omises ou seulement partiellement spécifiées: "cauchy", "gamma", "logistic", "binomial négatif" (paramétré par mu et taille), "t" et "weibull ". Notez que ces valeurs de départ peuvent ne pas être assez bonnes si l'ajustement est médiocre: en particulier, elles ne résistent pas aux valeurs aberrantes à moins que la distribution ajustée ne soit à longue queue.

Sergey Bushmanov
la source
1
Les deux réponses (Flom et Bushmanov) sont utiles. Je choisis celui-ci, car il rend plus explicite qu'avec les bonnes valeurs initiales et les contraintes, l'optimisation «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")

entrez la description de l'image ici

n

Essayons quelques simulations:

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

L'affichage de l'estimation est très instable (en regardant l'histogramme, une partie importante des valeurs estimées est à la limite supérieure donnée pour optimiser de 200).

Répéter avec un échantillon plus grand:

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

ce qui est beaucoup mieux, mais la moyenne est toujours bien supérieure à la valeur réelle de 2,5.

N'oubliez pas qu'il s'agit d'une version simplifiée du problème réel où les paramètres de localisation et d'échelle doivent également être estimés.

tν

kjetil b halvorsen
la source
5
Votre conclusion selon laquelle les problèmes d'estimation de df pourraient réellement aller à l'encontre de la raison de choisir une distribution t en premier lieu (c'est-à-dire la robustesse) est une réflexion qui suscite la réflexion.
user12719
1
(+1) "Unbounded above" n'est pas une mauvaise réponse et pourrait bien être utile à certaines fins lorsqu'il est associé à une estimation d'intervalle. L'important n'est pas d'utiliser aveuglément les informations de Fisher observées pour former des intervalles de confiance de Wald.
Scortchi - Réintégrer Monica
8

Dans l'aide pour fitdistr se trouve cet exemple:

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

indiquant que vous avez juste besoin d'une valeur pour df. Mais cela suppose une normalisation.

Pour plus de contrôle, ils montrent également

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))

où les paramètres seraient m = moyenne, s = écart type, df = degrés de liberté

Peter Flom - Réintégrer Monica
la source
1
Je suppose que je suis confus au sujet des paramètres d'une distribution en T. Y a-t-il 2 paramètres (moyenne, df) ou 3 (moyenne, écart type, df)? Je me demandais si l'on pouvait adapter le paramètre «df».
user12719
1
@ user12719 La distribution de Student-t a trois paramètres: emplacement, échelle et degrés de liberté. Ils ne sont pas appelés moyenne, écart-type et df car la moyenne et la variance de cette distribution dépendent des trois paramètres. De plus, ils n'existent pas dans certains cas. Peter Flom corrige le df mais cela peut également être considéré comme un paramètre inconnu.
1
@PeterFlom Dans le cas de la distribution de Cauchy, il est explicite que m et s sont l'emplacement et l'échelle. Je suis d'accord que les notations m et s suggèrent qu'elles représentent respectivement la moyenne et l'écart type. Mais cela peut être simplement une simplification de \muet \sigmaaussi. +1 il y a longtemps, au fait.
1
@PeterFlom Cette citation du fichier d'aide de R implique-t-elle que df est toujours 9 pour la distribution des étudiants? Ne pensez-vous pas que df devrait également être estimé? En fait, l'absence de dfest la cause de l'erreur, et la bonne réponse devrait fournir une recette pour la trouver.
Sergey Bushmanov
1
@PeterFlom BTW, si vous lisez le fichier d'aide quelques lignes au-dessus de votre citation, vous découvrirez pourquoi df=9est bon dans leur exemple et non pertinent ici.
Sergey Bushmanov, le