Please see PDF version
NIACSYNIA as a Tool for Statisticians
J. Michael Steele, Princeton University
ABSTRACI this report is an artifact of using the typesetting features
Three concrete examples are used to show how sym of MACSYMA ( check out the flag setting typeset:true).
bolic: computation can be of value in statistical research. (cl) difF(sin(x)/xx,3~,
Enough of the syntax of MACSYMA is discussed to pro~
vide an introduction to actual usage. The questions con
sidered concern the relationship of conditional cumulants 3 sin(x 6 sin(x cos(x ) + 6 cos(x (d 1)
Of homogeneous statistics, a method for proving limit X X X
theorems by studying zeros of polynomials, and the
Inendre transformation of some probability densities. (c2) taylor(cos(x)^ 3,xAQ
One question which is constantly considered is * When
will MACSYMA seriously help 3x 2 + 7x 4 _ 61xl + (d2)
1. Introduction 2 8
MACSYMA and other symbolic computational tools
are now widely available. Still, the statistical commun (c3) solve(2*x^ 2+x+2Y
ity has not felt the full impact of these new tools and
many people may not yet realize that the conduct of '115i + 1 J i~i 1
research in probability and statistics will be forever [X ,X 1 WM
changed by their presence. The purpose of this report is 4 4
to illustrate how even naive use of MACSYMA can be of
substantial benefit in research activity which is far from (cS) quit(Y,
being a priori computational.
The plan of the report is to first give a gentile intro
duction to the structure and syntax of MACSYMA, give This microsession is enough to show some of the
a few hints for using MACSYMA comfortably In a power in symbolic computation. One does not need
[INIX environment, and then give three examples. The MACSYMA to calculate the third derivative of sin(x)/x,
are not concocted but are just selected from the or the Taylor series Of COS(X )21 to four terms, but it is in
rny own experience in using MACSYMA as one of the fact easier and more reliable than doing it by hand.
tools of daytoday research. Finally, brief comments are Definite and indefinite integrals are also available with a
given con(erning the ways one can customize a comparably natural synt,m The use of solve illustrates a
\4;\('SYMA environment to facilitate the types of com slight variation on natural syntax but
putdtJons common in probability and statistics. solve(2*x 2+x+3)=0 works just as well.
Before going to a real example, one should learn
2. Getting Started some quick syntax. First "*,* terminates a command (or c
The way to get started using MACSYMA is to fire ) line. Also, [a,bA is a list containing a, b, and C. There
are assignment operators *." which defines functions
up MACSYMA, open An Introduction to MACSYMA
by Joel Moses and the MACSYMA Group of Synibloics~ and for more typical assignment.
Ind1984), and start playing. This is fun and effective. In 3. Some Practical Hints
!stark contrast, you can go to MACSYMA with a tough There are some frustrations which most newcomers
problem which you have not succeeded in doing by
hand, try to use MACSYMA for the first time and end to MACSYMA seem to meet and which are easily over
up meeting only frustration and likely failure. come once 'you know the trick*. First, the beginner
should avoid the N4ACSYMA editor. The best way to use
In teaching MACSYMA to my friends, the hardest MACS~ Is to work interactively only while explor
thing to get across always seems to be that MACS~ is
!rig and playing. Once a serious problem is to be studied
just a tool. A proper attitude Is created by thinking of the MACSYMA commands should be put in a file, say
MACSYMA as being like integration by parts. Surely
compfile, and executed in MACS~ with the com
integration by parts is a magnificent tool. It shines in mand batchCcomphie)
vital areas like the calculus of variations. It gives birth to The benefits of this process are (1) If an error is
the formulas of Gauss, Green, and Stokes. It is a weapon
made you can stop MACS~ with ^ Z and then just
4 first enoill in many parts of asymptoties. But. for all edit compfile width your favorite editor, (2) The code
of thw virlues, integration by parts will not solve dil you create becomes an asset which can be used in
prnblexiis, and neither w ill MAC.SYMA. The first task of
difrerent contexts after the traditional pattern of
the beginner is to shitn the inevitable disappointment of "modify rather than writew, and (3) The code can be
naive P ~*~(tations and to 4rticulate those areas of per conveniently quoted In papers, shared with others, etc.
sonal mk.erest where MACS1i'MA has a. realistic chance of
Please do not be mislead. There are many other
tools in MACSYMA which allow one to communicate
Fhewope of MA(SYMA and its interactive nature with files and with structures within MACSYMA. but.
are casv to see by looking over a brief but complete
for a while, these seem best left to vertuosi.
UNIX script of a trial session. In the first computations
we review, the syntax should be obvious if one bears in Two more practical hints ? Even if you do not wish
mind that the clines are entered by the user and the d to do ANY symbolic computation MACSYNIA can be a
lines are MACSYNIA's responses. One small point is that ('y~nd for Its FORTRAN and eqn facilities. in a nut
^ is used to denote exponentiation. The symbols W11 Itell 1.()R'rRAN(exi)ression) produces the FOR1RAN
W2), etc. will appear on the left side on the terminal just (ofle it) compute the specified expression. 1 fin cl this a
like the clabels. Their appearance on the right side in v(r~ pleasant way to avoid coding mistake., In the same
way the typeset flag of MACSYNIA permits, you to get 'he
t e.qn code of the resulting dline. This is much more pleasant and precise than writing your own equations with eqn. The equations of this report were all done in this way.
4. Simplicity 'rhough Diligence
One way in which MACSYNIA can be of aid in research is simply by removing the difficulties of looking at the problem from many sides. Often one has an idea which can go unexplored because the calculations look too cumbersome. This inhibition emerges for two reasons. First, it takes conviction to go trough a tedious calculation which may not lead to a useful result and second, the notion of usefulness is often linked to simplicity.
This block can be busted by MACSYMA in a number of cases. One from my experience concerns conditional cumulants of homogeneous statistics. The reason on can expect MACSYMA to be efFective here is that the tasks required are just tedious, not ingenious. At difFerentiation, substitution, and rational simplification by expansion MACSYMA is as surefooted as a mountain goat.
5. Cumulant Correction Formulas
We will assume that f satisfies the homogeneity condition (1.1) and set
Z', f(X1X2,,M (2.1)
where the Xi are independent and uniform on [0,11d. Further, we denote by Z, a random variable with the same definition as Z,, except that one of the Xi's is constrained to be on "forward" part of the boundary of [0, 1 ]d; that is, one of the Xi's is constrained to be in the Set B _ [0, Ild [0, 1)d.
If Y is any random variable the cumulant generating function of Y is defined by
7_Ki(Y) L log EW') (12)
i =1 i !
and the numbers Kl(Y) are the cumulants (or semiinvariants) of Y. Naturally, Kl(Y) is the mean Of Y, KZ is the variance, and these are the most important of the cumulants. The third and forth cumulants ( the skewness and the kurtosis) are of interest principally because after the variance all of the cumulants of a normal distribution are zero. The cumulants of Z, do not have to tend to zero for it to be asymptotically normal, but if the third and forth cumulants tend to zero it is a strong in(ii(.dlion of asymptotic normality.
One of our goals in pursuing the simulation method stimlied here (and in particular with bothering about higlier cumulants) is to obtain a compulationally ,.\1)~;dieni rnethod of assessing the possibility of asymptot ic noi inality in the case of computationally important random variables like the length of the minimal spanning trm In Rainey (1982) the asymptotic normality of the MST was proved under the hypothesis of an assumplion which seems very difficult to verify (but which is still quite plausible). The asymptotic normality of the length of the MST remains an important open problem. In fact, there is no subadditive Euclidean functional for which normality has been roved yet it is quite likely to hold for the whole class (` c.f. Beardwood, Halton, and flammersley (1959) and Steele (1981)).
To begin the conditioning argument, we let
S = minis: Xi E [O,Sld) thus, S is the size of the side of the smallest cube contain
!rig the data and the origin. We note that with probability one there is exactly one of the Xi on the forward part of the boundary of [O,Sld. A key observation is that we then have that the distribution of Z,, equals the distribution of SPZ,, where 5 and Z,, are generated independently.
Now we can just calculate.
exp K, (Z" ) 1, T ~ E(EW ~S)) (2.3)
Next we expand both sides of the above using
g 1 2
exp F a, t + a, + (a2 + a2)!,'
i =1 7! Z
+ (a 3 + 3a2a l + a3)iii +
(a 4 + 6a2a 2 + 4aa, + 3a 2 tAl
1 1 2 + a,) ~! +
and then complete taking the expectation by using the elementary fact that E (SO') dn (ot + dn )1. On simplifying the algebra, one finds the following tidy expressions for the conditional cumulants (To save space we will write Ki for Ki (Z, ) :
KI(Z',,) = K10 + L) (2.4a)
K2(Z', K2(1 + 2P) K 2 p (2.4b)
nd n 2 d
K3(7.'n) = K3(1 +
and finally we have
3p 6K 1K2p2 2K3 3
+ 1 p 2.4c)
nd nTdT illTI71
2 3 4 4
+ 24K1K2P 6KIP
2) P 2
12(K1K3 + KZ F 2.4d) n rd~
The method we have used is a natural one, but there are two surprises. The first of these is that the homogeneity hypothesis is powerful enough to give such simple results. The second is that it turns out to be incredibly neater to give an expression for the conditional variances in terms of the unconditional ones rather than viceversa.
The idea of using a generating function approach to obtain conditional cumulant formulas is inherent in the definition of cumulants. It has been used often, and in particular the method is applied in Brillinger(1969) to obtain a general formula relating cumulant to conditional cumulants. Interestingly, that an attempt to obtain the identities (2.4ad) directly from Brillinger's general formula is almost certainly doomed to failure because of the very much greater algebraic complexity of the expressions for the unconditional cumulants in terms of the conditional ones.
rhere is a second method for obtaining results relating conditional and inconditional cumulants which uses Mobius function of the lattice of partitions. Speed(I983) introduced this method and showed its attractiveness in a number of examples including an intriguing derivation
of Brillinger's formula. The Moblus function methoa For example, the following famous transform relation is
shows great promise; but, at the present level 1 of complex the key to Holders inequality,
ity, series inversion used here still seems to be the tool of
choice. XY < _XP + Y~
Some comments on the formulas (2.4ad) are in p q
order before going on to the simulation method. In the provided 1+ 1. In the same way that one can use
first place, even the modest equation (2.4b) is able to pro p q
vide some interesting lower bounds on the variance of Holder's inequality to determine the the dual space of
the pth power integrable functions, one can start off
certain subadditive Euclidean functionals. We will illus with any Legendre transform pair.
trate with the traveling salesman problem (TSP) in
dimension 2 , but comparable results will obviously hold A statistical problem of interest is the determination
for the minimal spanning tree, minimal triangulation, of all of those functions h (g,o) which have a finite
and other functionals and other dimensions. integral when multiplied by the normal density
The only fact we require about the TSP is that if (b(x ;g, a) and integrated with respect to g and a. To
(X lly 2, '. ' X, ) is the length of the shortest tour stud this MACSYMA was used to calculate the (bivari
through the points Xi 1 Q ~i then K1 is asymptotic to ate ~ Legendre transformation of (t with respect to
1 A and a. The result takes a page to express. Still, this
~n _f as n tends to infinity (Beardwood, Halton, Ham worm of an expression is probably the best available tool
for determining which priors on p and a will guarantee
mersley (1959)). Since this f is homogeneous of order
1 a posterior distribution which is absolutely continuous
p=l, equation (2.4b) implies that K2 is fX ). This modest with respect to Lebesgue measure.
bound is of interest because it is known that K2 IS For a discussion of the classical Legendre transfor
[x)unded for all n ( Steele (1981) ), but no better lower mation ( as it applies to Hamiltonian systems ) as well as
bound than the one just give is known. Naturally, it is one version of a MACSYMA Legendre transformation
expected for the TSP in d=2 that K2 converges to a posi one should consult Rand (1984).
tive constant as n goes to infinity. Mixed Computations
6. Simplicity Though Diligence The symbolic computational capacities Of
MACSYMA are augmented with numerical and graphi
One way in which MACSYMA can be of aid in cal capabilities. This addition to the workbench permits
research is simply by removing the difficulties of looking
the joint investigation of "formulas" and the exploration
at the problem from many sides. Often one has an idea of their quantitative features.
which can go unexplored because the calculations look
too cumbersome. This inhibition emerges for two reasons. One way this plays out to the investigators advan
First, it takes conviction to go trough a tedious calcula tage is that it is possible to gain insight from much larger
tion which may not lead to a useful result, and second, examples than would normally be possible. A useful
example of this from my own experience come for the
the notion of usefulness is often linked to simplicity. theory of random trees.
This block can be busted by MACSYMA in a The detailed example is given in Steele(1985), but
number of cases. One from my experience concerns con the role of MACSYMA is easy to review without those
ditional cumulants of homogeneous statistics. The reason details. The generating function for the number of leaves
on can expect MACSYMA to be efective here is that the of a random tree was derived in Renyi(I959) under the
tasks required are just tedious, not ingenious. At model of randomness which corresponds to the uniform
differentiation, substitution, and rational simplification
distribution on the Pruffer codes ( These codes provide a
by expansion MACSYMA is as surefooted as a mountain natural onetoone correspondence between labeled
unrooted trees and n2 tuples of the integers 1 <, k ~n
One of the properties of mathematics ( which also For the purpose of modeling, the uniform model on
persists in statistics) is a strong preference for simple
trees is quite limited , eg. one can not escape the conse
~ers. Part of the process ot getting to understand the
power of symbolic computation rest in comming to terms quence of there being approximately leaves on the
with complex results. What makes this novel is that one great majority of the trees. flow much nicer for model
typically does not have to deal with complex results ing purposes if we had a family of tree valued random
because once a calculation gets past a certain complexity
variables whose expected number of trees could be
threshold it is given up for lost. In a computational chosen by the modeler ? There is such a family and a
environment one can get extremely complex results and method for generating choices from that family given in
one faces the new challenge of using such results intelli
One virtue of Renyi's distribution on random trees
is that Reny! was able to prove a central limit theorem
the number of leaves of a tree form his family ( as n
7. The Legendre Experience tends to infinity). To establish the corresponding result
in the exponential family of random trees it appeared to
This feature of the changing role of complexity is be almost necessary to show that the generating func
well illustrated by the enjoyable practice of computing tions given by Renyi had all real zeros. Naturally that
Inendre transformations. For any function f(x) the opinion takes considerable motivation. Using MACSYNIA
Lnendre transformation may be defined by it was possible to examine polynomials up to n=10 were
g (y max (xy f (x D it was discovered that for large n there were nontrivial
complex roots. The plan of analysis was saved when it
The importance of this transformation in analysis was found that by adding a stochasticly negligible term
comes In part from the relaxed inequality
to Renyi's polynomial the roots could be forced back to
ry < f (X) + g (y) the real line. This eventually lead to a proof of the cen
tral limit theorem for the number of leaves of a tree Gladd, N. T. (1984) Symbolic Manipulation Techniques
from the exponential family. This is particularly satisfy for Plasma Kinetic Theory Derivations Journal of Com
ing because it plays back into the original motivation by Putational Physics, 56, 175202.
making the mean number of leaves an essential parame
ter. Hochbaum, D. and Steele J. M. (19821 Steinhaus's
Q. Concluding Remarks geometric location problem for random samples In the
1 This introduction to MACSYMA has been casual plane, Adv. App4 Prob, 14, 5667.
and discursive. The goal has been to motivate others to Karp, R. M. (1977). Probabalistic analysis of partitioning
begin productive interaction with the tools which are algorithms for the traveling salesman problem in the
now widely available. The hints which have been given plane. Math. Oper. Res. 2 202224.
are those that seem to be responsible and welcoming.
There is a great deal of benefit to the statistical commun MACS~ Group of Symbolics, IncPrincipal writer
ity which can be gained from symbolic computation and Joel Moses(I984), An Introduction to MACSYMA Sym
it is inevitable that current elforts will be made to look bolics. Inc. Cambridge Mass.
primative by those which will follow very shortly. It Is
especially intersting to imagine how our collective con
ception of what constitutes a useful answer will migrate Mathlab Group (1983) MACSYMA Reference Manual,
in the next few years. version Ten, Vol. I & 11, Laboratory for Computer Science
ACKNOWLEDGEMENT MIT, Cambridge~ Man
This research was supported in part by NSF grant
DMS8414069. Results reproduced here were obtained Ramey. D. (1982) A Nonparametric Test of Bimodality
in part by using MACSYMA, a large symbolic manipula with Application to Cluster Analysis. Department of
tion program developed at the MIT laboratory for com StatisticsYale University, New Haven Conn.
puter Science and supported from 1975 to 1983 by Rand, R. FL (1984) Computer Algebra in Applied
NASA under grant NSG 1323. by ONR under grant
N0001477C0641. by DOE under grant ET78C02 Mathematics. an Introduction to MACSYMA Pitman
4687,by the USAF under F4%2079C020. and since Advanced Publishing Program, Boston.
1982 by Symbolics, Inc. of Cambridge, Mass. MACSYMA
is a trademark of Symbolics, Inc. Key Words. Reny~ A. and Sulanke, RX19641 Uber die konvexe Hulle
MACSYMA, symbolic computation, cumulants, Legendre von n zufallig gewahIten Punkten. IL Z Wahrscheinli
transformation, central limit theorem. AMS Subject chtkeitsteorie, 3 138147
Classifi cations (19801 Major 68C20; Minor 62EX
Author's Address: Program in Engineering Statistics . Renyi, A~1956). Some Remarks on the Theory of Trees.
Fngineering Quadrangle, Princeton University, Princeton Pub. Math. Inst. Hung. Acad. Sci. 4, 7385.
N. J. 08544 Speed, T. (1983). Cumulants and Partition Lattices. Aust.
J. Statist. 25, 378388
REFERENCES Steele, J. NI. (1980). Subaditive Euclidean Functionals
and nonlinear growth in geometric probability. Ann.
Beardwood, J. Halton, H. J. and Harnmersley, 1M. (1959). Probability. 2 7585.
The shortest path through many ponts. Proc. Camb.
Phil. Soc. 55, 299327. Steele, J. M. (1985). Exponential Family of Trees. Techni
cal Report, Program in Engineering Statistics and
Brillinger, D.R. (1969). The calculation of cumulants via Management Science, Princeton University.
conditioning. Ann. Inst. Statist. Math., 21, 375390.