Implémentation de Dirichlet cdf?

8

J'ai besoin de calculer le CDF de Dirichlet , mais je ne trouve que les implémentations du PDF .

Connaissez-vous une bibliothèque (de préférence dans R) l'implémentant?

Ricky Robinson
la source
1
Pas directement au courant. Mais il peut y avoir quelque chose qui peut être fait. Que devez-vous en faire?
Glen_b -Reinstate Monica
2
Je dois prendre la complémentarité du CDF et le considérer comme ma valeur p.
Ricky Robinson
2
Hmm. Alors ouais, si tu as besoin1P(X1x1,X2x2,...,Xkxk), vous avez en quelque sorte besoin du cdf. L'idée de simulation de Zen est certainement une façon de le faire (et plus le nombre de dimensions est élevé, mieux il commence à paraître), mais si vous le faites, utilisez l'un des packages avec les implémentations intégrées de rdirichlet. S'il ne s'agit que de 3 ou peut-être de 4 (le dernier composant, bien sûr, étant redondant), cela peut valoir la peine d'essayer la quadrature numérique.
Glen_b -Reinstate Monica

Réponses:

9

N'oubliez pas que si Yi sont indépendants Gamma(ai,b), pour i=1,,k, puis

(X1,,Xk)=(Y1j=1kYj,,Ykj=1kYj)Dirichlet(a1,,ak).

La preuve se trouve à la page 594 du livre de Luc Devroye .

Par conséquent, une possibilité consiste à calculer une approximation de Monte Carlo de

FX1,,Xk(t1,,tk)=P{X1t1,,Xktk},
en commençant par les gammas. Dans R, essayez ceci:
pdirichlet <- function(a, t) {
    N <- 10000
    rdirichlet <- function(a) { y <- rgamma(length(a), a, 1); y / sum(y) }
    x <- replicate(N, rdirichlet(a), simplify = FALSE)
    sum(sapply(x, function(x) prod(x <= t))) / N
}

Je n'ai pas vérifié le code. Utilisez-le soigneusement. Si vous trouvez des bogues, veuillez nous en informer.

Zen
la source
3
Il y a déjà une fonction rdirichlet vectorisée dans R - en fait plusieurs d'entre elles (dans gtools, MCMCpacket dirmultpar exemple).
Glen_b -Reinstate Monica
@Zen J'ai a <- c(6, 20,2)comment obtenir le cdf Drichelt? est la matrice t 2 par 2?
score324
J'ai utilisé le code ci-dessus, mais cela génère une erreur.
score324
@ score324 Un vecteur Dirichlet a généralement chaque élément dans [0,1] et leur somme étant 1, donc (6,20,2)est en dehors du support
Henry
2

Une bibliothèque? Mathematica l'a. Voici le code d'un exemple de tracé d'un CDF Dirichlet de la documentation:

Plot3D[CDF[DirichletDistribution[{1, 3, 2}], {x, y}], {x, 0, 1}, {y, 0, 1}]
Mike Z.
la source
1
Comment obtenir l'expression pour CDF [DirichletDistribution [{1, 3, 2}]]?
AIB