rim - R’s interface to Maxima
Development version of the CRAN package. rim provides an
interface to the powerful and fairly complete maxima computer algebra
system
This repository is a fork from the original RMaxima
created by Kseniia Shumelchyk and Hans W. Borcher shumelchyk/rmaxima,
which is currently not maintained.
Requirements
- Maxima (>= 5.42.1), this package has been tested with version between 5.42.1 and 5.46.0, but may also work with older versios
Installation
To install the current release version from CRAN:
install.packages("rim")If you want to install the latest version the easiest way is to
install the R package remotes first and install the package
from this github repo:
install.packages("remotes")
remotes::install_github(repo = "rcst/rim")Usage
If you want to learn how to use Maxima you can check out a number of ressources from it’s project site.
RMarkdown
This section demonstrates the use of rim via it’s knitr
engine. Alternatively, you can run each line of the Maxima code chunks
via maxima.get().
This page has been generated using rim’s knitr engine
itself.
Note that you can set the output format,
e.g. maxima.options(engine.format = "latex") to get
LaTeX-style expression rendering. When using pandoc to render the
resulting file, MathJax is loaded or embedded automatically, when the
ouput is a HTML document and so the Maxima output expressions are
rendered nicely. Alternatively, expression outputs can be printed in
MathML directly using
maxima.options(engine.format = "mathml").
When attaching the package and Maxima is installed already, then it
will be set up and started in the background automatically, when a
command is send. In addition, if you have installed the R-package
knitr (which get’s installed by default with this package
as well), it will register maxima as a knitr engine. If
this worked, maxima code can be used and rendered in RMarkdown
documents.
# library(rim)
devtools::load_all()## ℹ Loading rim
## Maxima successfully registered as knitr engine!maxima.options(engine.format = "latex",
engine.label = TRUE,
inline.format = "latex",
inline.label = FALSE)For example, to generate this page you can download it’s source
index.Rmd from here
and then use rmarkdown::render to generate it.
download.file(url = "https://raw.githubusercontent.com/rcst/rim/master/docs/index.Rmd",
destfile = "index.Rmd")
rmarkdown::render("index.Rmd")Now we can enter Maxima expression inside code chunks. Note that we
need to end each line by ; or $ (suppressing
output). For example we can type the following code chunk into our
RMarkdown document file:
```{maxima}
L: sqrt(1 - 1/R^2);
assume(R > 0)$
'integrate(x, x, 0, L) = integrate(x, x, 0, L);
```
In the above code chunk we define a variable L (that
depends on another variable R). We tell Maximas database to
assume R being larger than zero and suppress any output
resulting from this command. The last line prints an unevaluated
integral on the left-hand side of the equal sign and on the right-hand
side evaluates that same definit integral. The tick quotation mark tells
Maxima not to evaluate the suceeding expression. This, when rendered
into a document will be printed as
(%i1) L: sqrt(1 - 1/R^2);\[\mathtt{(\textit{%o}_{1})}\quad \sqrt{1-\frac{1}{R^2}}\]
(%i2) assume(R > 0)$
(%i3) 'integrate(x, x, 0, L) = integrate(x, x, 0, L);\[\mathtt{(\textit{%o}_{3})}\quad \int_{0}^{\sqrt{1-\frac{1}{R^2}}}{x\;dx}=\frac{R^2-1}{2\,R^2}\]
The resulting output (if not suppressed) is printed after each line of code, and thus apprears as several. Input and output reference labels that are assign by Maxima are printed into he code chunk. Those can be used to refer to previous commands in other code chunks. For example
(%i4) sqrt(rhs(%o3));\[\mathtt{(\textit{%o}_{4})}\quad \frac{\sqrt{R^2-1}}{\sqrt{2}\,R}\]
takes the right-hand side of the result of the equation from input
label %i3, which is assigned to output label
%o3 and computes the square-root of it.
Of course exercising on such simple computations is of little benefit. The real benefit comes from more complicated expression and the effort that we would need to put if we wanted to typeset the result, such as this
(%i5) integrate(1 / (1 + x^4), x);\[\mathtt{(\textit{%o}_{5})}\quad \frac{\log \left(x^2+\sqrt{2}\,x+1\right)}{2^{\frac{5}{2}}}-\frac{\log \left(x^2-\sqrt{2}\,x+1\right)}{2^{\frac{5}{2}}}+\frac{\arctan \left(\frac{2\,x+\sqrt{2}}{\sqrt{2}}\right)}{2^{\frac{3}{2}}}+\frac{\arctan \left(\frac{2\,x-\sqrt{2}}{\sqrt{2}}\right)}{2^{\frac{3}{2}}}\]
pandoc automatically renders the LaTeX output format
from Maxima by including MathJax JavaScript script. In general,
pandoc takes care of how mathematical equations delimited
by $$ are rendered.
However, we can also change Maxima’s output format to MathML, which works if the output document is a HTML document.
maxima.options(engine.format = "mathml")
maxima.options()## $format
## [1] "linear"
##
## $engine.format
## [1] "mathml"
##
## $inline.format
## [1] "inline"
##
## $label
## [1] TRUE
##
## $engine.label
## [1] TRUE
##
## $inline.label
## [1] FALSE(%i6) sqrt(3/4);(%i7) f(x) := e^(x^2)$
(%i8) diff(f(x), x);(%i9) %;Notice, that we can use the symbol % to refer to the
last expression and that this works across code chunks.
You can also replay other previous expressions by referring to their
output labels, as already demonstrated above. Whether or not reference
labels are printed via the knitr engine has no influence on
their existence in Maxima internally.
rim also supports printing Maxima output as inline. This
works by putting an r inline code chunk and using
maxima.inline(). Example:
`r maxima.inline("log(%o1);")` results in \(\frac{\log
\left(1-\frac{1}{R^2}\right)}{2}\) .
maxima.options(engine.format = "latex")Comments can be added to Maxima code as text between /*
and */ (example taken from the Maxima
manual).
(%i11) /* aa is a variable of interest */ aa : 1234;\[\mathtt{(\textit{%o}_{11})}\quad 1234\]
(%i12) bb : aa^2; /* Value of bb depends on aa */ \[\mathtt{(\textit{%o}_{12})}\quad 1522756\]
(%i13) /* User-defined infix operator */ infix ("Q");\[\mathtt{(\textit{%o}_{13})}\quad \mbox{ Q }\]
(%i14) /* Parses same as a b c, not abc */ a/* foo */Q/* bar */c;\[\mathtt{(\textit{%o}_{14})}\quad \textit{aQc}\]
(%i15) /* Comments /* can be nested /* to any depth */ */ */ 1 + xyz;\[\mathtt{(\textit{%o}_{15})}\quad \textit{xyz}+1\]
Chunk Options
The following chunk options are currently available:
echo: Should the code chunk be printed (default isTRUE). Ifeval=FALSE, individual code lines are printed without without input labels, since these are assigned by Maxima upon evaluationeval: Should the code chunk be evaluated by Maxima (default isTRUE)include: whether code chunk and results should be printed (default isTRUE). IfFALSE, code is still evaluatedoutput.var: A character string for a variable that caputures the parsed output expression(s) (i.e., translated into R-code), see example below.
NOTE: The results chunk option is
explicitly ignored from chunk header and implicitly controlled via
maxima.options(engine.format = ...).
Capturing parsed Maxima output
One of the above Maxima code chunk was actually set with option
output.var = "foo"
```{maxima output-var, output.var = "foo"}
integrate(1 / (1 + x^4), x);
```
which captured the Maxima result(s) as parsed R-code (unevaluated expressions) into a named list. Those expressions can subsequently be evaluated in R e.g., to evaluate an expression in the context of a data.frame
df <- data.frame(x = 1:10, y = letters[1:10])
df$z <- eval(foo[["o5"]], envir = df)
df## x y z
## 1 1 a 0.866973
## 2 2 b 1.070128
## 3 3 c 1.098440
## 4 4 d 1.105521
## 5 5 e 1.108056
## 6 6 f 1.109178
## 7 7 g 1.109749
## 8 8 h 1.110070
## 9 9 i 1.110264
## 10 10 j 1.110387foo## $o5
## (((2L^(-3L/2L)) * atan(((2L^(-1L/2L)) * ((-1L * (2L^(1L/2L))) +
## (2L * x))))) + ((2L^(-3L/2L)) * atan(((2L^(-1L/2L)) * ((2L^(1L/2L)) +
## (2L * x))))) + (-1L * (2L^(-5L/2L)) * log((1L + (-1L * (2L^(1L/2L)) *
## x) + (x^2L)))) + ((2L^(-5L/2L)) * log((1L + ((2L^(1L/2L)) *
## x) + (x^2L)))))In this way, Maxima can assist us with symolic computations in our statistical modelling in R.
Demos
Plots
(%i16) r: (exp(cos(t))-2*cos(4*t)-sin(t/12)^5)$
(%i17) plot2d([parametric, r*sin(t), r*cos(t), [t,-8*%pi,8*%pi]]);
plot2d()
(%i18) plot3d(log (x^2*y^2), [x, -2, 2], [y, -2, 2],[grid, 29, 29],
[palette, [gradient, red, orange, yellow, green]],
color_bar, [xtics, 1], [ytics, 1], [ztics, 4],
[color_bar_tics, 4]);
plot3d()
(%i19) example1:
gr3d (title = "Controlling color range",
enhanced3d = true,
color = green,
cbrange = [-3,10],
explicit(x^2+y^2, x,-2,2,y,-2,2)) $
(%i20) example2:
gr3d (title = "Playing with tics in colorbox",
enhanced3d = true,
color = green,
cbtics = {["High",10],["Medium",05],["Low",0]},
cbrange = [0, 10],
explicit(x^2+y^2, x,-2,2,y,-2,2))$
(%i21) example3:
gr3d (title = "Logarithmic scale to colors",
enhanced3d = true,
color = green,
logcb = true,
logz = true,
palette = [-15,24,-9],
explicit(exp(x^2-y^2), x,-2,2,y,-2,2))$
(%i22) draw(
dimensions = [500,1500],
example1, example2, example3);
draw()
(%i23) draw2d(
dimensions = [1000, 1000],
proportional_axes = xy,
fill_color = sea_green,
color = aquamarine,
line_width = 6,
ellipse(7,6,2,3,0,360));
draw2d()
(%i24) draw3d(
dimensions = [1000, 1000],
surface_hide = true,
axis_3d = false,
proportional_axes = xyz,
color = blue,
cylindrical(z,z,-2,2,a,0,2*%pi),
color = brown,
cylindrical(3,z,-2,2,az,0,%pi),
color = green,
cylindrical(sqrt(25-z^2),z,-5,5,a,0,%pi));
draw3d()
Distributions
This is basically a reproduced demo for Sympy from Bryan
Zhang’s Blog.
First we define some helper functions:
Normal Distribution
(%i30) normal(x) :=
(2*%pi*sigma^2)^(-1/2) *
exp(-(x-mu)^2/(2*sigma^2));\[\mathtt{(\textit{%o}_{30})}\quad \textit{normal}\left(x\right):=\left(2\,\pi\,\sigma^2\right)^{\frac{-1}{2}}\,\exp \left(\frac{-\left(x-\mu\right)^2}{2\,\sigma^2}\right)\]
(%i31) assume(sigma > 0)$
(%i32) area(normal(x));\[\mathtt{(\textit{%o}_{32})}\quad 1\]
(%i33) mean(normal(x));\[\mathtt{(\textit{%o}_{33})}\quad \mu\]
(%i34) variance(normal(x));\[\mathtt{(\textit{%o}_{34})}\quad \frac{2^{\frac{3}{2}}\,\sqrt{\pi}\,\sigma^3+2^{\frac{3}{2}}\,\sqrt{\pi}\,\mu^2\,\sigma}{2^{\frac{3}{2}}\,\sqrt{\pi}\,\sigma}-\mu^2\]
(%i35) mgf(normal(x));\[\mathtt{(\textit{%o}_{35})}\quad e^{\frac{\sigma^2\,t^2+2\,\mu\,t}{2}}\]
Pareto Distribution
(%i36) pareto(x) := (alpha * (x[m])^alpha) / x^(alpha+1);\[\mathtt{(\textit{%o}_{36})}\quad \textit{pareto}\left(x\right):=\frac{\alpha\,x_{m}^{\alpha}}{x^{\alpha+1}}\]
Laplace Distribution
(%i39) laplace(x) := (2*b)^-1 * exp(-abs(x - mu)/b);\[\mathtt{(\textit{%o}_{39})}\quad \textit{laplace}\left(x\right):=\left(2\,b\right)^ {- 1 }\,\exp \left(\frac{-\left| x-\mu\right| }{b}\right)\]
(%i40) load("abs_integrate")$
(%i41) assume(b > 0)$
(%i42) area(laplace(x));\[\mathtt{(\textit{%o}_{42})}\quad 1\]
(%i43) mean(laplace(x));\[\mathtt{(\textit{%o}_{43})}\quad \mu\]
(%i44) variance(laplace(x));\[\mathtt{(\textit{%o}_{44})}\quad \frac{2\,b\,\mu^2+4\,b^3}{2\,b}-\mu^2\]
Exponential Distribution
(%i45) expo(x) := unit_step(x) * lambda * exp(-lambda * x);\[\mathtt{(\textit{%o}_{45})}\quad \textit{expo}\left(x\right):=\textit{unit\_step}\left(x\right)\,\lambda\,\exp \left(\left(-\lambda\right)\,x\right)\]
(%i46) assume(lambda > 0)$
(%i47) area(expo(x));\[\mathtt{(\textit{%o}_{47})}\quad 1\]
(%i48) mean(expo(x));\[\mathtt{(\textit{%o}_{48})}\quad \frac{1}{\lambda}\]
(%i49) variance(expo(x));\[\mathtt{(\textit{%o}_{49})}\quad \frac{1}{\lambda^2}\]
Matrices
(%i50) m: matrix([0, 1, a], [1, 0, 1], [1, 1, 0]);\[\mathtt{(\textit{%o}_{50})}\quad \begin{pmatrix}0 & 1 & a \\ 1 & 0 & 1 \\ 1 & 1 & 0 \\ \end{pmatrix}\]
(%i51) transpose(m);\[\mathtt{(\textit{%o}_{51})}\quad \begin{pmatrix}0 & 1 & 1 \\ 1 & 0 & 1 \\ a & 1 & 0 \\ \end{pmatrix}\]
(%i52) determinant(m);\[\mathtt{(\textit{%o}_{52})}\quad a+1\]
(%i53) f: invert(m), detout;\[\mathtt{(\textit{%o}_{53})}\quad \frac{\begin{pmatrix}-1 & a & 1 \\ 1 & -a & a \\ 1 & 1 & -1 \\ \end{pmatrix}}{a+1}\]
(%i54) m . f;\[\mathtt{(\textit{%o}_{54})}\quad \begin{pmatrix}0 & 1 & a \\ 1 & 0 & 1 \\ 1 & 1 & 0 \\ \end{pmatrix}\cdot \left(\frac{\begin{pmatrix}-1 & a & 1 \\ 1 & -a & a \\ 1 & 1 & -1 \\ \end{pmatrix}}{a+1}\right)\]
(%i55) expand(%);\[\mathtt{(\textit{%o}_{55})}\quad \begin{pmatrix}\frac{a}{a+1}+\frac{1}{a+1} & 0 & 0 \\ 0 & \frac{a}{a+1}+\frac{1}{a+1} & 0 \\ 0 & 0 & \frac{a}{a+1}+\frac{1}{a+1} \\ \end{pmatrix}\]
(%i56) factor(%);\[\mathtt{(\textit{%o}_{56})}\quad \begin{pmatrix}1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \\ \end{pmatrix}\]
If-then-else
(%i57) x: 1234;\[\mathtt{(\textit{%o}_{57})}\quad 1234\]
(%i58) y: 2345;\[\mathtt{(\textit{%o}_{58})}\quad 2345\]
(%i59) if x > y
then x
else y;\[\mathtt{(\textit{%o}_{59})}\quad 2345\]