I've been trying to evaluate sage for integration and expression
manipulation in the context of probabilistic models. I understand it
uses Maxima for this. After posting on the sage mailing list I was
pointed here for further help/discussion.
My problem is that Maxima can integrate a beta distribution but not a
Dirichlet distribution.
Here's my original post on the sage list:
http://tinyurl.com/3yb7wld
As discussed I'd be willing to try to implement some aspect of this,
assuming I'm capable enough and given a few pointers on how to start.
Here's my sage code:
a, b, c = var('a b c')
assume(a > 0)
assume(b > 0)
assume(c > 0)
#
# Beta
#
x = var('x')
beta_dist = x**(a-1) * (1 - x)**(b-1)
c = integral(beta_dist, x, 0, 1)
#
# Dirichlet
#
x_1, x_2, x_3 = var('x_1, x_2, x_3')
dirichlet_dist = x_1**(a-1) * x_2**(b-1) * (1-x_2-x_1)**(c-1)
e = integral(dirichlet_dist, x_1, 0, 1)
f = integral(e, x_2, 0, 1)
Then in my sage session I attach the code and get:
----------------------------------------------------------------------
| Sage Version 4.4.2, Release Date: 2010-05-19 |
| Type notebook() for the GUI, and license() for information. |
----------------------------------------------------------------------
sage: attach('hdpm.py')
sage: c
beta(a, b)
sage: f
integrate(x_2^(b - 1)*integrate((-x_1 - x_2 + 1)^(beta(a, b) - 1)*x_1^(a
- 1), x_1, 0, 1), x_2, 0, 1)
Presumably it left the integrals alone because it couldn't do them. Is
there any way to help it along?