File: stanmath.html

package info (click to toggle)
r-cran-stanheaders 2.32.10-1
links: PTS, VCS
area: main
in suites: forky, sid, trixie
size: 19,652 kB
sloc: cpp: 125,776; ansic: 73,706; f90: 30,959; ml: 243; sh: 15; makefile: 5
file content (1342 lines) | stat: -rw-r--r-- 127,825 bytes
<!DOCTYPE html>

<html>

<head>

<meta charset="utf-8" />
<meta name="generator" content="pandoc" />
<meta http-equiv="X-UA-Compatible" content="IE=EDGE" />

<meta name="viewport" content="width=device-width, initial-scale=1" />

<meta name="author" content="Stan Development Team" />

<meta name="date" content="2024-07-15" />

<title>Using the Stan Math C++ Library</title>


<style type="text/css">
  code{white-space: pre-wrap;}
  span.smallcaps{font-variant: small-caps;}
  span.underline{text-decoration: underline;}
  div.column{display: inline-block; vertical-align: top; width: 50%;}
  div.hanging-indent{margin-left: 1.5em; text-indent: -1.5em;}
  ul.task-list{list-style: none;}
    </style>



<style type="text/css">
  code {
    white-space: pre;
  }
  .sourceCode {
    overflow: visible;
  }
</style>
<style type="text/css" data-origin="pandoc">
a.sourceLine { display: inline-block; line-height: 1.25; }
a.sourceLine { pointer-events: none; color: inherit; text-decoration: inherit; }
a.sourceLine:empty { height: 1.2em; }
.sourceCode { overflow: visible; }
code.sourceCode { white-space: pre; position: relative; }
div.sourceCode { margin: 1em 0; }
pre.sourceCode { margin: 0; }
@media screen {
div.sourceCode { overflow: auto; }
}
@media print {
code.sourceCode { white-space: pre-wrap; }
a.sourceLine { text-indent: -1em; padding-left: 1em; }
}
pre.numberSource a.sourceLine
  { position: relative; left: -4em; }
pre.numberSource a.sourceLine::before
  { content: attr(title);
    position: relative; left: -1em; text-align: right; vertical-align: baseline;
    border: none; pointer-events: all; display: inline-block;
    -webkit-touch-callout: none; -webkit-user-select: none;
    -khtml-user-select: none; -moz-user-select: none;
    -ms-user-select: none; user-select: none;
    padding: 0 4px; width: 4em;
    color: #aaaaaa;
  }
pre.numberSource { margin-left: 3em; border-left: 1px solid #aaaaaa;  padding-left: 4px; }
div.sourceCode
  {  }
@media screen {
a.sourceLine::before { text-decoration: underline; }
}
code span.al { color: #ff0000; font-weight: bold; } /* Alert */
code span.an { color: #60a0b0; font-weight: bold; font-style: italic; } /* Annotation */
code span.at { color: #7d9029; } /* Attribute */
code span.bn { color: #40a070; } /* BaseN */
code span.bu { } /* BuiltIn */
code span.cf { color: #007020; font-weight: bold; } /* ControlFlow */
code span.ch { color: #4070a0; } /* Char */
code span.cn { color: #880000; } /* Constant */
code span.co { color: #60a0b0; font-style: italic; } /* Comment */
code span.cv { color: #60a0b0; font-weight: bold; font-style: italic; } /* CommentVar */
code span.do { color: #ba2121; font-style: italic; } /* Documentation */
code span.dt { color: #902000; } /* DataType */
code span.dv { color: #40a070; } /* DecVal */
code span.er { color: #ff0000; font-weight: bold; } /* Error */
code span.ex { } /* Extension */
code span.fl { color: #40a070; } /* Float */
code span.fu { color: #06287e; } /* Function */
code span.im { } /* Import */
code span.in { color: #60a0b0; font-weight: bold; font-style: italic; } /* Information */
code span.kw { color: #007020; font-weight: bold; } /* Keyword */
code span.op { color: #666666; } /* Operator */
code span.ot { color: #007020; } /* Other */
code span.pp { color: #bc7a00; } /* Preprocessor */
code span.sc { color: #4070a0; } /* SpecialChar */
code span.ss { color: #bb6688; } /* SpecialString */
code span.st { color: #4070a0; } /* String */
code span.va { color: #19177c; } /* Variable */
code span.vs { color: #4070a0; } /* VerbatimString */
code span.wa { color: #60a0b0; font-weight: bold; font-style: italic; } /* Warning */

</style>
<script>
// apply pandoc div.sourceCode style to pre.sourceCode instead
(function() {
  var sheets = document.styleSheets;
  for (var i = 0; i < sheets.length; i++) {
    if (sheets[i].ownerNode.dataset["origin"] !== "pandoc") continue;
    try { var rules = sheets[i].cssRules; } catch (e) { continue; }
    var j = 0;
    while (j < rules.length) {
      var rule = rules[j];
      // check if there is a div.sourceCode rule
      if (rule.type !== rule.STYLE_RULE || rule.selectorText !== "div.sourceCode") {
        j++;
        continue;
      }
      var style = rule.style.cssText;
      // check if color or background-color is set
      if (rule.style.color === '' && rule.style.backgroundColor === '') {
        j++;
        continue;
      }
      // replace div.sourceCode by a pre.sourceCode rule
      sheets[i].deleteRule(j);
      sheets[i].insertRule('pre.sourceCode{' + style + '}', j);
    }
  }
})();
</script>




<style type="text/css">body {
background-color: #fff;
margin: 1em auto;
max-width: 700px;
overflow: visible;
padding-left: 2em;
padding-right: 2em;
font-family: "Open Sans", "Helvetica Neue", Helvetica, Arial, sans-serif;
font-size: 14px;
line-height: 1.35;
}
#TOC {
clear: both;
margin: 0 0 10px 10px;
padding: 4px;
width: 400px;
border: 1px solid #CCCCCC;
border-radius: 5px;
background-color: #f6f6f6;
font-size: 13px;
line-height: 1.3;
}
#TOC .toctitle {
font-weight: bold;
font-size: 15px;
margin-left: 5px;
}
#TOC ul {
padding-left: 40px;
margin-left: -1.5em;
margin-top: 5px;
margin-bottom: 5px;
}
#TOC ul ul {
margin-left: -2em;
}
#TOC li {
line-height: 16px;
}
table {
margin: 1em auto;
border-width: 1px;
border-color: #DDDDDD;
border-style: outset;
border-collapse: collapse;
}
table th {
border-width: 2px;
padding: 5px;
border-style: inset;
}
table td {
border-width: 1px;
border-style: inset;
line-height: 18px;
padding: 5px 5px;
}
table, table th, table td {
border-left-style: none;
border-right-style: none;
}
table thead, table tr.even {
background-color: #f7f7f7;
}
p {
margin: 0.5em 0;
}
blockquote {
background-color: #f6f6f6;
padding: 0.25em 0.75em;
}
hr {
border-style: solid;
border: none;
border-top: 1px solid #777;
margin: 28px 0;
}
dl {
margin-left: 0;
}
dl dd {
margin-bottom: 13px;
margin-left: 13px;
}
dl dt {
font-weight: bold;
}
ul {
margin-top: 0;
}
ul li {
list-style: circle outside;
}
ul ul {
margin-bottom: 0;
}
pre, code {
background-color: #f7f7f7;
border-radius: 3px;
color: #333;
white-space: pre-wrap; 
}
pre {
border-radius: 3px;
margin: 5px 0px 10px 0px;
padding: 10px;
}
pre:not([class]) {
background-color: #f7f7f7;
}
code {
font-family: Consolas, Monaco, 'Courier New', monospace;
font-size: 85%;
}
p > code, li > code {
padding: 2px 0px;
}
div.figure {
text-align: center;
}
img {
background-color: #FFFFFF;
padding: 2px;
border: 1px solid #DDDDDD;
border-radius: 3px;
border: 1px solid #CCCCCC;
margin: 0 5px;
}
h1 {
margin-top: 0;
font-size: 35px;
line-height: 40px;
}
h2 {
border-bottom: 4px solid #f7f7f7;
padding-top: 10px;
padding-bottom: 2px;
font-size: 145%;
}
h3 {
border-bottom: 2px solid #f7f7f7;
padding-top: 10px;
font-size: 120%;
}
h4 {
border-bottom: 1px solid #f7f7f7;
margin-left: 8px;
font-size: 105%;
}
h5, h6 {
border-bottom: 1px solid #ccc;
font-size: 105%;
}
a {
color: #0033dd;
text-decoration: none;
}
a:hover {
color: #6666ff; }
a:visited {
color: #800080; }
a:visited:hover {
color: #BB00BB; }
a[href^="http:"] {
text-decoration: underline; }
a[href^="https:"] {
text-decoration: underline; }

code > span.kw { color: #555; font-weight: bold; } 
code > span.dt { color: #902000; } 
code > span.dv { color: #40a070; } 
code > span.bn { color: #d14; } 
code > span.fl { color: #d14; } 
code > span.ch { color: #d14; } 
code > span.st { color: #d14; } 
code > span.co { color: #888888; font-style: italic; } 
code > span.ot { color: #007020; } 
code > span.al { color: #ff0000; font-weight: bold; } 
code > span.fu { color: #900; font-weight: bold; } 
code > span.er { color: #a61717; background-color: #e3d2d2; } 
</style>




</head>

<body>




<h1 class="title toc-ignore">Using the Stan Math C++ Library</h1>
<h4 class="author">Stan Development Team</h4>
<h4 class="date">2024-07-15</h4>



<div id="using-the-stanheaders-package-from-other-r-packages" class="section level1">
<h1>Using the <strong>StanHeaders</strong> Package from Other R Packages</h1>
<p>The <strong>StanHeaders</strong> package contains no R functions. To use the Stan Math Library in other packages, it is often sufficient to specify</p>
<pre><code>LinkingTo: StanHeaders (&gt;= 2.26.0), RcppParallel (&gt;= 5.0.1)</code></pre>
<p>in the DESCRIPTION file of another package and put something like</p>
<pre><code>CXX_STD = CXX17
PKG_CXXFLAGS = $(shell &quot;$(R_HOME)/bin$(R_ARCH_BIN)/Rscript&quot; -e &quot;RcppParallel::CxxFlags()&quot;) \
               $(shell &quot;$(R_HOME)/bin$(R_ARCH_BIN)/Rscript&quot; -e &quot;StanHeaders:::CxxFlags()&quot;)
PKG_LIBS = $(shell &quot;$(R_HOME)/bin$(R_ARCH_BIN)/Rscript&quot; -e &quot;RcppParallel::RcppParallelLibs()&quot;) \
           $(shell &quot;$(R_HOME)/bin$(R_ARCH_BIN)/Rscript&quot; -e &quot;StanHeaders:::LdFlags()&quot;)</code></pre>
<p>in the src/Makevars and src/Makevars.win files and put <code>GNU make</code> in the <code>SystemRequirements:</code> field of the package’s DESCRIPTION file. If, in addition, the other package needs to utilize the MCMC, optimization, variational inference, or parsing facilities of the Stan Library, then it is also necessary to include the <code>src</code> directory of <strong>StanHeaders</strong> in the other package’s <code>PKG_CXXFLAGS</code> in the src/Makevars and src/Makevars.win files with something like</p>
<pre><code>STANHEADERS_SRC = $(shell &quot;$(R_HOME)/bin$(R_ARCH_BIN)/Rscript&quot; -e &quot;message()&quot; \
  -e &quot;cat(system.file('include', 'src', package = 'StanHeaders', mustWork = TRUE))&quot; \
  -e &quot;message()&quot; | grep &quot;StanHeaders&quot;)
PKG_CXXFLAGS += -I&quot;$(STANHEADERS_SRC)&quot;</code></pre>
</div>
<div id="calling-functions-in-the-stanheaders-package-from-r" class="section level1">
<h1>Calling functions in the <strong>StanHeaders</strong> Package from R</h1>
<p>The only exposed R function in the in the <strong>StanHeaders</strong> package is <code>stanFunction</code>, which can be used to call most functions in the Stan Math Library.</p>
<div class="sourceCode" id="cb4"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb4-1" title="1"><span class="kw">example</span>(stanFunction, <span class="dt">package =</span> <span class="st">&quot;StanHeaders&quot;</span>, <span class="dt">run.dontrun =</span> <span class="ot">TRUE</span>)</a>
<a class="sourceLine" id="cb4-2" title="2"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-3" title="3"><span class="co">#&gt; stnFnc&gt;   files &lt;- dir(system.file(&quot;include&quot;, &quot;stan&quot;, &quot;math&quot;, &quot;prim&quot;,</span></a>
<a class="sourceLine" id="cb4-4" title="4"><span class="co">#&gt; stnFnc+                            package = &quot;StanHeaders&quot;), </span></a>
<a class="sourceLine" id="cb4-5" title="5"><span class="co">#&gt; stnFnc+                pattern = &quot;hpp$&quot;, recursive = TRUE)</span></a>
<a class="sourceLine" id="cb4-6" title="6"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-7" title="7"><span class="co">#&gt; stnFnc&gt;   functions &lt;- sub(&quot;\\.hpp$&quot;, &quot;&quot;, </span></a>
<a class="sourceLine" id="cb4-8" title="8"><span class="co">#&gt; stnFnc+                    sort(unique(basename(files[dirname(files) != &quot;.&quot;]))))</span></a>
<a class="sourceLine" id="cb4-9" title="9"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-10" title="10"><span class="co">#&gt; stnFnc&gt;   length(functions) # you could call most of these Stan functions</span></a>
<a class="sourceLine" id="cb4-11" title="11"><span class="co">#&gt; [1] 955</span></a>
<a class="sourceLine" id="cb4-12" title="12"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-13" title="13"><span class="co">#&gt; stnFnc&gt;     log(sum(exp(exp(1)), exp(pi))) # true value</span></a>
<a class="sourceLine" id="cb4-14" title="14"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-15" title="15"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-16" title="16"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;log_sum_exp&quot;, x = exp(1), y = pi)</span></a>
<a class="sourceLine" id="cb4-17" title="17"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-18" title="18"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-19" title="19"><span class="co">#&gt; stnFnc&gt;     args(log_sum_exp) # now exists in .GlobalEnv</span></a>
<a class="sourceLine" id="cb4-20" title="20"><span class="co">#&gt; function (x, y) </span></a>
<a class="sourceLine" id="cb4-21" title="21"><span class="co">#&gt; NULL</span></a>
<a class="sourceLine" id="cb4-22" title="22"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-23" title="23"><span class="co">#&gt; stnFnc&gt;     log_sum_exp(x = pi, y = exp(1))</span></a>
<a class="sourceLine" id="cb4-24" title="24"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-25" title="25"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-26" title="26"><span class="co">#&gt; stnFnc&gt;     # but log_sum_exp() was not defined for a vector or matrix</span></a>
<a class="sourceLine" id="cb4-27" title="27"><span class="co">#&gt; stnFnc&gt;     x &lt;- c(exp(1), pi)</span></a>
<a class="sourceLine" id="cb4-28" title="28"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-29" title="29"><span class="co">#&gt; stnFnc&gt;     try(log_sum_exp(x))</span></a>
<a class="sourceLine" id="cb4-30" title="30"><span class="co">#&gt; Error in log_sum_exp(x) : argument &quot;y&quot; is missing, with no default</span></a>
<a class="sourceLine" id="cb4-31" title="31"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-32" title="32"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;log_sum_exp&quot;, x = x) # now it is</span></a>
<a class="sourceLine" id="cb4-33" title="33"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-34" title="34"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-35" title="35"><span class="co">#&gt; stnFnc&gt;     # log_sum_exp() is now also defined for a matrix</span></a>
<a class="sourceLine" id="cb4-36" title="36"><span class="co">#&gt; stnFnc&gt;     log_sum_exp(as.matrix(x))</span></a>
<a class="sourceLine" id="cb4-37" title="37"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-38" title="38"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-39" title="39"><span class="co">#&gt; stnFnc&gt;     log_sum_exp(t(as.matrix(x)))</span></a>
<a class="sourceLine" id="cb4-40" title="40"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-41" title="41"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-42" title="42"><span class="co">#&gt; stnFnc&gt;     log_sum_exp(rbind(x, x))</span></a>
<a class="sourceLine" id="cb4-43" title="43"><span class="co">#&gt; [1] 4.338465</span></a>
<a class="sourceLine" id="cb4-44" title="44"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-45" title="45"><span class="co">#&gt; stnFnc&gt;     # but log_sum_exp() was not defined for a list</span></a>
<a class="sourceLine" id="cb4-46" title="46"><span class="co">#&gt; stnFnc&gt;     try(log_sum_exp(as.list(x)))</span></a>
<a class="sourceLine" id="cb4-47" title="47"><span class="co">#&gt; Error in eval(ei, envir) : </span></a>
<a class="sourceLine" id="cb4-48" title="48"><span class="co">#&gt;   Not compatible with requested type: [type=list; target=double].</span></a>
<a class="sourceLine" id="cb4-49" title="49"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-50" title="50"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;log_sum_exp&quot;, x = as.list(x)) # now it is</span></a>
<a class="sourceLine" id="cb4-51" title="51"><span class="co">#&gt; [1] 3.645318</span></a>
<a class="sourceLine" id="cb4-52" title="52"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-53" title="53"><span class="co">#&gt; stnFnc&gt;     # in rare cases, passing a nested list is needed</span></a>
<a class="sourceLine" id="cb4-54" title="54"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;dims&quot;, x = list(list(1:3)))</span></a>
<a class="sourceLine" id="cb4-55" title="55"><span class="co">#&gt; [1] 1 1 3</span></a>
<a class="sourceLine" id="cb4-56" title="56"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-57" title="57"><span class="co">#&gt; stnFnc&gt;     # functions of complex arguments work</span></a>
<a class="sourceLine" id="cb4-58" title="58"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;eigenvalues&quot;, # different ordering than base:eigen()</span></a>
<a class="sourceLine" id="cb4-59" title="59"><span class="co">#&gt; stnFnc+                  x = matrix(complex(real = 1:9, imaginary = pi),</span></a>
<a class="sourceLine" id="cb4-60" title="60"><span class="co">#&gt; stnFnc+                             nrow = 3, ncol = 3))</span></a>
<a class="sourceLine" id="cb4-61" title="61"><span class="co">#&gt; [1] -8.179880e-17+1.622085e-18i -8.612657e-01+4.854073e-01i  1.586127e+01+8.939371e+00i</span></a>
<a class="sourceLine" id="cb4-62" title="62"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-63" title="63"><span class="co">#&gt; stnFnc&gt;     # nullary functions work but are not that interesting</span></a>
<a class="sourceLine" id="cb4-64" title="64"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;negative_infinity&quot;)</span></a>
<a class="sourceLine" id="cb4-65" title="65"><span class="co">#&gt; [1] -Inf</span></a>
<a class="sourceLine" id="cb4-66" title="66"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-67" title="67"><span class="co">#&gt; stnFnc&gt;     # PRNG functions work by adding a seed argument</span></a>
<a class="sourceLine" id="cb4-68" title="68"><span class="co">#&gt; stnFnc&gt;     stanFunction(&quot;lkj_corr_rng&quot;, K = 3L, eta = 1)</span></a>
<a class="sourceLine" id="cb4-69" title="69"><span class="co">#&gt;            [,1]      [,2]       [,3]</span></a>
<a class="sourceLine" id="cb4-70" title="70"><span class="co">#&gt; [1,] 1.00000000 0.5654780 0.03415387</span></a>
<a class="sourceLine" id="cb4-71" title="71"><span class="co">#&gt; [2,] 0.56547796 1.0000000 0.24813894</span></a>
<a class="sourceLine" id="cb4-72" title="72"><span class="co">#&gt; [3,] 0.03415387 0.2481389 1.00000000</span></a>
<a class="sourceLine" id="cb4-73" title="73"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb4-74" title="74"><span class="co">#&gt; stnFnc&gt;     args(lkj_corr_rng) # has a seed argument</span></a>
<a class="sourceLine" id="cb4-75" title="75"><span class="co">#&gt; function (K, eta, random_seed = sample.int(.Machine$integer.max, </span></a>
<a class="sourceLine" id="cb4-76" title="76"><span class="co">#&gt;     size = 1L)) </span></a>
<a class="sourceLine" id="cb4-77" title="77"><span class="co">#&gt; NULL</span></a></code></pre></div>
<style type="text/css">
.scroll-100 {
  max-height: 100px;
  overflow-y: auto;
  background-color: inherit;
}
</style>
<p>The <code>functions</code> object defined in this example lists the many Stan functions that could be called (if all their arguments are numeric, see <code>help(stanFunction, package = &quot;StanHeaders&quot;)</code> for details)</p>
<pre><code>#&gt;        [,1]                                [,2]                                        
#&gt;   [1,] &quot;Eigen&quot;                             &quot;LDLT_factor&quot;                               
#&gt;   [2,] &quot;Phi&quot;                               &quot;Phi_approx&quot;                                
#&gt;   [3,] &quot;StdVectorBuilder&quot;                  &quot;VectorBuilder&quot;                             
#&gt;   [4,] &quot;VectorBuilderHelper&quot;               &quot;abs&quot;                                       
#&gt;   [5,] &quot;accumulator&quot;                       &quot;acos&quot;                                      
#&gt;   [6,] &quot;acosh&quot;                             &quot;ad_promotable&quot;                             
#&gt;   [7,] &quot;add&quot;                               &quot;add_diag&quot;                                  
#&gt;   [8,] &quot;algebra_solver_adapter&quot;            &quot;all&quot;                                       
#&gt;   [9,] &quot;any&quot;                               &quot;append_array&quot;                              
#&gt;  [10,] &quot;append_col&quot;                        &quot;append_return_type&quot;                        
#&gt;  [11,] &quot;append_row&quot;                        &quot;apply&quot;                                     
#&gt;  [12,] &quot;apply_scalar_binary&quot;               &quot;apply_scalar_ternary&quot;                      
#&gt;  [13,] &quot;apply_scalar_unary&quot;                &quot;apply_vector_unary&quot;                        
#&gt;  [14,] &quot;arg&quot;                               &quot;array_builder&quot;                             
#&gt;  [15,] &quot;as_array_or_scalar&quot;                &quot;as_bool&quot;                                   
#&gt;  [16,] &quot;as_column_vector_or_scalar&quot;        &quot;as_value_array_or_scalar&quot;                  
#&gt;  [17,] &quot;as_value_column_array_or_scalar&quot;   &quot;as_value_column_vector_or_scalar&quot;          
#&gt;  [18,] &quot;asin&quot;                              &quot;asinh&quot;                                     
#&gt;  [19,] &quot;assign&quot;                            &quot;atan&quot;                                      
#&gt;  [20,] &quot;atan2&quot;                             &quot;atanh&quot;                                     
#&gt;  [21,] &quot;autocorrelation&quot;                   &quot;autocovariance&quot;                            
#&gt;  [22,] &quot;base_type&quot;                         &quot;bernoulli_ccdf_log&quot;                        
#&gt;  [23,] &quot;bernoulli_cdf&quot;                     &quot;bernoulli_cdf_log&quot;                         
#&gt;  [24,] &quot;bernoulli_lccdf&quot;                   &quot;bernoulli_lcdf&quot;                            
#&gt;  [25,] &quot;bernoulli_log&quot;                     &quot;bernoulli_logit_glm_log&quot;                   
#&gt;  [26,] &quot;bernoulli_logit_glm_lpmf&quot;          &quot;bernoulli_logit_glm_rng&quot;                   
#&gt;  [27,] &quot;bernoulli_logit_log&quot;               &quot;bernoulli_logit_lpmf&quot;                      
#&gt;  [28,] &quot;bernoulli_logit_rng&quot;               &quot;bernoulli_lpmf&quot;                            
#&gt;  [29,] &quot;bernoulli_rng&quot;                     &quot;bessel_first_kind&quot;                         
#&gt;  [30,] &quot;bessel_second_kind&quot;                &quot;beta&quot;                                      
#&gt;  [31,] &quot;beta_binomial_ccdf_log&quot;            &quot;beta_binomial_cdf&quot;                         
#&gt;  [32,] &quot;beta_binomial_cdf_log&quot;             &quot;beta_binomial_lccdf&quot;                       
#&gt;  [33,] &quot;beta_binomial_lcdf&quot;                &quot;beta_binomial_log&quot;                         
#&gt;  [34,] &quot;beta_binomial_lpmf&quot;                &quot;beta_binomial_rng&quot;                         
#&gt;  [35,] &quot;beta_ccdf_log&quot;                     &quot;beta_cdf&quot;                                  
#&gt;  [36,] &quot;beta_cdf_log&quot;                      &quot;beta_lccdf&quot;                                
#&gt;  [37,] &quot;beta_lcdf&quot;                         &quot;beta_log&quot;                                  
#&gt;  [38,] &quot;beta_lpdf&quot;                         &quot;beta_proportion_ccdf_log&quot;                  
#&gt;  [39,] &quot;beta_proportion_cdf_log&quot;           &quot;beta_proportion_lccdf&quot;                     
#&gt;  [40,] &quot;beta_proportion_lcdf&quot;              &quot;beta_proportion_log&quot;                       
#&gt;  [41,] &quot;beta_proportion_lpdf&quot;              &quot;beta_proportion_rng&quot;                       
#&gt;  [42,] &quot;beta_rng&quot;                          &quot;binary_log_loss&quot;                           
#&gt;  [43,] &quot;binomial_ccdf_log&quot;                 &quot;binomial_cdf&quot;                              
#&gt;  [44,] &quot;binomial_cdf_log&quot;                  &quot;binomial_coefficient_log&quot;                  
#&gt;  [45,] &quot;binomial_lccdf&quot;                    &quot;binomial_lcdf&quot;                             
#&gt;  [46,] &quot;binomial_log&quot;                      &quot;binomial_logit_log&quot;                        
#&gt;  [47,] &quot;binomial_logit_lpmf&quot;               &quot;binomial_lpmf&quot;                             
#&gt;  [48,] &quot;binomial_rng&quot;                      &quot;block&quot;                                     
#&gt;  [49,] &quot;bool_constant&quot;                     &quot;boost_policy&quot;                              
#&gt;  [50,] &quot;broadcast_array&quot;                   &quot;categorical_log&quot;                           
#&gt;  [51,] &quot;categorical_logit_glm_lpmf&quot;        &quot;categorical_logit_log&quot;                     
#&gt;  [52,] &quot;categorical_logit_lpmf&quot;            &quot;categorical_logit_rng&quot;                     
#&gt;  [53,] &quot;categorical_lpmf&quot;                  &quot;categorical_rng&quot;                           
#&gt;  [54,] &quot;cauchy_ccdf_log&quot;                   &quot;cauchy_cdf&quot;                                
#&gt;  [55,] &quot;cauchy_cdf_log&quot;                    &quot;cauchy_lccdf&quot;                              
#&gt;  [56,] &quot;cauchy_lcdf&quot;                       &quot;cauchy_log&quot;                                
#&gt;  [57,] &quot;cauchy_lpdf&quot;                       &quot;cauchy_rng&quot;                                
#&gt;  [58,] &quot;cbrt&quot;                              &quot;ceil&quot;                                      
#&gt;  [59,] &quot;check_2F1_converges&quot;               &quot;check_3F2_converges&quot;                       
#&gt;  [60,] &quot;check_bounded&quot;                     &quot;check_cholesky_factor&quot;                     
#&gt;  [61,] &quot;check_cholesky_factor_corr&quot;        &quot;check_column_index&quot;                        
#&gt;  [62,] &quot;check_consistent_size&quot;             &quot;check_consistent_sizes&quot;                    
#&gt;  [63,] &quot;check_consistent_sizes_mvt&quot;        &quot;check_corr_matrix&quot;                         
#&gt;  [64,] &quot;check_cov_matrix&quot;                  &quot;check_finite&quot;                              
#&gt;  [65,] &quot;check_flag_sundials&quot;               &quot;check_greater&quot;                             
#&gt;  [66,] &quot;check_greater_or_equal&quot;            &quot;check_ldlt_factor&quot;                         
#&gt;  [67,] &quot;check_less&quot;                        &quot;check_less_or_equal&quot;                       
#&gt;  [68,] &quot;check_lower_triangular&quot;            &quot;check_matching_dims&quot;                       
#&gt;  [69,] &quot;check_matching_sizes&quot;              &quot;check_multiplicable&quot;                       
#&gt;  [70,] &quot;check_nonnegative&quot;                 &quot;check_nonzero_size&quot;                        
#&gt;  [71,] &quot;check_not_nan&quot;                     &quot;check_ordered&quot;                             
#&gt;  [72,] &quot;check_pos_definite&quot;                &quot;check_pos_semidefinite&quot;                    
#&gt;  [73,] &quot;check_positive&quot;                    &quot;check_positive_finite&quot;                     
#&gt;  [74,] &quot;check_positive_ordered&quot;            &quot;check_range&quot;                               
#&gt;  [75,] &quot;check_row_index&quot;                   &quot;check_simplex&quot;                             
#&gt;  [76,] &quot;check_size_match&quot;                  &quot;check_sorted&quot;                              
#&gt;  [77,] &quot;check_square&quot;                      &quot;check_std_vector_index&quot;                    
#&gt;  [78,] &quot;check_symmetric&quot;                   &quot;check_unit_vector&quot;                         
#&gt;  [79,] &quot;check_vector&quot;                      &quot;check_vector_index&quot;                        
#&gt;  [80,] &quot;chi_square_ccdf_log&quot;               &quot;chi_square_cdf&quot;                            
#&gt;  [81,] &quot;chi_square_cdf_log&quot;                &quot;chi_square_lccdf&quot;                          
#&gt;  [82,] &quot;chi_square_lcdf&quot;                   &quot;chi_square_log&quot;                            
#&gt;  [83,] &quot;chi_square_lpdf&quot;                   &quot;chi_square_rng&quot;                            
#&gt;  [84,] &quot;child_type&quot;                        &quot;chol2inv&quot;                                  
#&gt;  [85,] &quot;cholesky_corr_constrain&quot;           &quot;cholesky_corr_free&quot;                        
#&gt;  [86,] &quot;cholesky_decompose&quot;                &quot;cholesky_factor_constrain&quot;                 
#&gt;  [87,] &quot;cholesky_factor_free&quot;              &quot;choose&quot;                                    
#&gt;  [88,] &quot;col&quot;                               &quot;cols&quot;                                      
#&gt;  [89,] &quot;columns_dot_product&quot;               &quot;columns_dot_self&quot;                          
#&gt;  [90,] &quot;compiler_attributes&quot;               &quot;complex_base&quot;                              
#&gt;  [91,] &quot;complex_schur_decompose&quot;           &quot;conj&quot;                                      
#&gt;  [92,] &quot;conjunction&quot;                       &quot;constants&quot;                                 
#&gt;  [93,] &quot;constraint_tolerance&quot;              &quot;contains_fvar&quot;                             
#&gt;  [94,] &quot;contains_std_vector&quot;               &quot;copysign&quot;                                  
#&gt;  [95,] &quot;corr_constrain&quot;                    &quot;corr_free&quot;                                 
#&gt;  [96,] &quot;corr_matrix_constrain&quot;             &quot;corr_matrix_free&quot;                          
#&gt;  [97,] &quot;cos&quot;                               &quot;cosh&quot;                                      
#&gt;  [98,] &quot;coupled_ode_system&quot;                &quot;cov_exp_quad&quot;                              
#&gt;  [99,] &quot;cov_matrix_constrain&quot;              &quot;cov_matrix_constrain_lkj&quot;                  
#&gt; [100,] &quot;cov_matrix_free&quot;                   &quot;cov_matrix_free_lkj&quot;                       
#&gt; [101,] &quot;crossprod&quot;                         &quot;csr_extract&quot;                               
#&gt; [102,] &quot;csr_extract_u&quot;                     &quot;csr_extract_v&quot;                             
#&gt; [103,] &quot;csr_extract_w&quot;                     &quot;csr_matrix_times_vector&quot;                   
#&gt; [104,] &quot;csr_to_dense_matrix&quot;               &quot;csr_u_to_z&quot;                                
#&gt; [105,] &quot;cumulative_sum&quot;                    &quot;determinant&quot;                               
#&gt; [106,] &quot;diag_matrix&quot;                       &quot;diag_post_multiply&quot;                        
#&gt; [107,] &quot;diag_pre_multiply&quot;                 &quot;diagonal&quot;                                  
#&gt; [108,] &quot;digamma&quot;                           &quot;dims&quot;                                      
#&gt; [109,] &quot;dirichlet_log&quot;                     &quot;dirichlet_lpdf&quot;                            
#&gt; [110,] &quot;dirichlet_lpmf&quot;                    &quot;dirichlet_rng&quot;                             
#&gt; [111,] &quot;discrete_range_ccdf_log&quot;           &quot;discrete_range_cdf&quot;                        
#&gt; [112,] &quot;discrete_range_cdf_log&quot;            &quot;discrete_range_lccdf&quot;                      
#&gt; [113,] &quot;discrete_range_lcdf&quot;               &quot;discrete_range_log&quot;                        
#&gt; [114,] &quot;discrete_range_lpmf&quot;               &quot;discrete_range_rng&quot;                        
#&gt; [115,] &quot;disjunction&quot;                       &quot;distance&quot;                                  
#&gt; [116,] &quot;divide&quot;                            &quot;divide_columns&quot;                            
#&gt; [117,] &quot;domain_error&quot;                      &quot;domain_error_vec&quot;                          
#&gt; [118,] &quot;dot&quot;                               &quot;dot_product&quot;                               
#&gt; [119,] &quot;dot_self&quot;                          &quot;double_exponential_ccdf_log&quot;               
#&gt; [120,] &quot;double_exponential_cdf&quot;            &quot;double_exponential_cdf_log&quot;                
#&gt; [121,] &quot;double_exponential_lccdf&quot;          &quot;double_exponential_lcdf&quot;                   
#&gt; [122,] &quot;double_exponential_log&quot;            &quot;double_exponential_lpdf&quot;                   
#&gt; [123,] &quot;double_exponential_rng&quot;            &quot;eigen_comparisons&quot;                         
#&gt; [124,] &quot;eigendecompose&quot;                    &quot;eigendecompose_sym&quot;                        
#&gt; [125,] &quot;eigenvalues&quot;                       &quot;eigenvalues_sym&quot;                           
#&gt; [126,] &quot;eigenvectors&quot;                      &quot;eigenvectors_sym&quot;                          
#&gt; [127,] &quot;elementwise_check&quot;                 &quot;elt_divide&quot;                                
#&gt; [128,] &quot;elt_multiply&quot;                      &quot;erf&quot;                                       
#&gt; [129,] &quot;erfc&quot;                              &quot;error_index&quot;                               
#&gt; [130,] &quot;eval&quot;                              &quot;exp&quot;                                       
#&gt; [131,] &quot;exp2&quot;                              &quot;exp_mod_normal_ccdf_log&quot;                   
#&gt; [132,] &quot;exp_mod_normal_cdf&quot;                &quot;exp_mod_normal_cdf_log&quot;                    
#&gt; [133,] &quot;exp_mod_normal_lccdf&quot;              &quot;exp_mod_normal_lcdf&quot;                       
#&gt; [134,] &quot;exp_mod_normal_log&quot;                &quot;exp_mod_normal_lpdf&quot;                       
#&gt; [135,] &quot;exp_mod_normal_rng&quot;                &quot;expm1&quot;                                     
#&gt; [136,] &quot;exponential_ccdf_log&quot;              &quot;exponential_cdf&quot;                           
#&gt; [137,] &quot;exponential_cdf_log&quot;               &quot;exponential_lccdf&quot;                         
#&gt; [138,] &quot;exponential_lcdf&quot;                  &quot;exponential_log&quot;                           
#&gt; [139,] &quot;exponential_lpdf&quot;                  &quot;exponential_rng&quot;                           
#&gt; [140,] &quot;fabs&quot;                              &quot;factor_U&quot;                                  
#&gt; [141,] &quot;factor_cov_matrix&quot;                 &quot;falling_factorial&quot;                         
#&gt; [142,] &quot;fdim&quot;                              &quot;fft&quot;                                       
#&gt; [143,] &quot;fill&quot;                              &quot;finite_diff_gradient&quot;                      
#&gt; [144,] &quot;finite_diff_gradient_auto&quot;         &quot;finite_diff_stepsize&quot;                      
#&gt; [145,] &quot;floor&quot;                             &quot;fma&quot;                                       
#&gt; [146,] &quot;fmax&quot;                              &quot;fmin&quot;                                      
#&gt; [147,] &quot;fmod&quot;                              &quot;for_each&quot;                                  
#&gt; [148,] &quot;forward_as&quot;                        &quot;frechet_ccdf_log&quot;                          
#&gt; [149,] &quot;frechet_cdf&quot;                       &quot;frechet_cdf_log&quot;                           
#&gt; [150,] &quot;frechet_lccdf&quot;                     &quot;frechet_lcdf&quot;                              
#&gt; [151,] &quot;frechet_log&quot;                       &quot;frechet_lpdf&quot;                              
#&gt; [152,] &quot;frechet_rng&quot;                       &quot;gamma_ccdf_log&quot;                            
#&gt; [153,] &quot;gamma_cdf&quot;                         &quot;gamma_cdf_log&quot;                             
#&gt; [154,] &quot;gamma_lccdf&quot;                       &quot;gamma_lcdf&quot;                                
#&gt; [155,] &quot;gamma_log&quot;                         &quot;gamma_lpdf&quot;                                
#&gt; [156,] &quot;gamma_p&quot;                           &quot;gamma_q&quot;                                   
#&gt; [157,] &quot;gamma_rng&quot;                         &quot;gaussian_dlm_obs_log&quot;                      
#&gt; [158,] &quot;gaussian_dlm_obs_lpdf&quot;             &quot;gaussian_dlm_obs_rng&quot;                      
#&gt; [159,] &quot;generalized_inverse&quot;               &quot;get&quot;                                       
#&gt; [160,] &quot;get_base1&quot;                         &quot;get_base1_lhs&quot;                             
#&gt; [161,] &quot;get_imag&quot;                          &quot;get_lp&quot;                                    
#&gt; [162,] &quot;get_real&quot;                          &quot;gp_dot_prod_cov&quot;                           
#&gt; [163,] &quot;gp_exp_quad_cov&quot;                   &quot;gp_exponential_cov&quot;                        
#&gt; [164,] &quot;gp_matern32_cov&quot;                   &quot;gp_matern52_cov&quot;                           
#&gt; [165,] &quot;gp_periodic_cov&quot;                   &quot;grad_2F1&quot;                                  
#&gt; [166,] &quot;grad_F32&quot;                          &quot;grad_inc_beta&quot;                             
#&gt; [167,] &quot;grad_pFq&quot;                          &quot;grad_reg_inc_beta&quot;                         
#&gt; [168,] &quot;grad_reg_inc_gamma&quot;                &quot;grad_reg_lower_inc_gamma&quot;                  
#&gt; [169,] &quot;gumbel_ccdf_log&quot;                   &quot;gumbel_cdf&quot;                                
#&gt; [170,] &quot;gumbel_cdf_log&quot;                    &quot;gumbel_lccdf&quot;                              
#&gt; [171,] &quot;gumbel_lcdf&quot;                       &quot;gumbel_log&quot;                                
#&gt; [172,] &quot;gumbel_lpdf&quot;                       &quot;gumbel_rng&quot;                                
#&gt; [173,] &quot;hcubature&quot;                         &quot;head&quot;                                      
#&gt; [174,] &quot;hmm_check&quot;                         &quot;hmm_hidden_state_prob&quot;                     
#&gt; [175,] &quot;hmm_latent_rng&quot;                    &quot;hmm_marginal&quot;                              
#&gt; [176,] &quot;holder&quot;                            &quot;hypergeometric_2F1&quot;                        
#&gt; [177,] &quot;hypergeometric_2F2&quot;                &quot;hypergeometric_3F2&quot;                        
#&gt; [178,] &quot;hypergeometric_log&quot;                &quot;hypergeometric_lpmf&quot;                       
#&gt; [179,] &quot;hypergeometric_pFq&quot;                &quot;hypergeometric_rng&quot;                        
#&gt; [180,] &quot;hypot&quot;                             &quot;i_times&quot;                                   
#&gt; [181,] &quot;identity_constrain&quot;                &quot;identity_free&quot;                             
#&gt; [182,] &quot;identity_matrix&quot;                   &quot;if_else&quot;                                   
#&gt; [183,] &quot;imag&quot;                              &quot;inc_beta&quot;                                  
#&gt; [184,] &quot;inc_beta_dda&quot;                      &quot;inc_beta_ddb&quot;                              
#&gt; [185,] &quot;inc_beta_ddz&quot;                      &quot;include_summand&quot;                           
#&gt; [186,] &quot;index_apply&quot;                       &quot;index_type&quot;                                
#&gt; [187,] &quot;init_threadpool_tbb&quot;               &quot;initialize&quot;                                
#&gt; [188,] &quot;initialize_fill&quot;                   &quot;int_step&quot;                                  
#&gt; [189,] &quot;integrate_1d&quot;                      &quot;integrate_1d_adapter&quot;                      
#&gt; [190,] &quot;integrate_ode_rk45&quot;                &quot;integrate_ode_std_vector_interface_adapter&quot;
#&gt; [191,] &quot;inv&quot;                               &quot;inv_Phi&quot;                                   
#&gt; [192,] &quot;inv_chi_square_ccdf_log&quot;           &quot;inv_chi_square_cdf&quot;                        
#&gt; [193,] &quot;inv_chi_square_cdf_log&quot;            &quot;inv_chi_square_lccdf&quot;                      
#&gt; [194,] &quot;inv_chi_square_lcdf&quot;               &quot;inv_chi_square_log&quot;                        
#&gt; [195,] &quot;inv_chi_square_lpdf&quot;               &quot;inv_chi_square_rng&quot;                        
#&gt; [196,] &quot;inv_cloglog&quot;                       &quot;inv_erfc&quot;                                  
#&gt; [197,] &quot;inv_gamma_ccdf_log&quot;                &quot;inv_gamma_cdf&quot;                             
#&gt; [198,] &quot;inv_gamma_cdf_log&quot;                 &quot;inv_gamma_lccdf&quot;                           
#&gt; [199,] &quot;inv_gamma_lcdf&quot;                    &quot;inv_gamma_log&quot;                             
#&gt; [200,] &quot;inv_gamma_lpdf&quot;                    &quot;inv_gamma_rng&quot;                             
#&gt; [201,] &quot;inv_inc_beta&quot;                      &quot;inv_logit&quot;                                 
#&gt; [202,] &quot;inv_sqrt&quot;                          &quot;inv_square&quot;                                
#&gt; [203,] &quot;inv_wishart_cholesky_lpdf&quot;         &quot;inv_wishart_cholesky_rng&quot;                  
#&gt; [204,] &quot;inv_wishart_log&quot;                   &quot;inv_wishart_lpdf&quot;                          
#&gt; [205,] &quot;inv_wishart_rng&quot;                   &quot;invalid_argument&quot;                          
#&gt; [206,] &quot;invalid_argument_vec&quot;              &quot;inverse&quot;                                   
#&gt; [207,] &quot;inverse_softmax&quot;                   &quot;inverse_spd&quot;                               
#&gt; [208,] &quot;is_any_nan&quot;                        &quot;is_arena_matrix&quot;                           
#&gt; [209,] &quot;is_autodiff&quot;                       &quot;is_base_pointer_convertible&quot;               
#&gt; [210,] &quot;is_cholesky_factor&quot;                &quot;is_cholesky_factor_corr&quot;                   
#&gt; [211,] &quot;is_column_index&quot;                   &quot;is_complex&quot;                                
#&gt; [212,] &quot;is_constant&quot;                       &quot;is_container&quot;                              
#&gt; [213,] &quot;is_container_or_var_matrix&quot;        &quot;is_corr_matrix&quot;                            
#&gt; [214,] &quot;is_dense_dynamic&quot;                  &quot;is_detected&quot;                               
#&gt; [215,] &quot;is_double_or_int&quot;                  &quot;is_eigen&quot;                                  
#&gt; [216,] &quot;is_eigen_dense_base&quot;               &quot;is_eigen_dense_dynamic&quot;                    
#&gt; [217,] &quot;is_eigen_matrix&quot;                   &quot;is_eigen_matrix_base&quot;                      
#&gt; [218,] &quot;is_eigen_sparse_base&quot;              &quot;is_fvar&quot;                                   
#&gt; [219,] &quot;is_inf&quot;                            &quot;is_integer&quot;                                
#&gt; [220,] &quot;is_kernel_expression&quot;              &quot;is_ldlt_factor&quot;                            
#&gt; [221,] &quot;is_less_or_equal&quot;                  &quot;is_lower_triangular&quot;                       
#&gt; [222,] &quot;is_mat_finite&quot;                     &quot;is_matching_dims&quot;                          
#&gt; [223,] &quot;is_matching_size&quot;                  &quot;is_matrix&quot;                                 
#&gt; [224,] &quot;is_matrix_cl&quot;                      &quot;is_nan&quot;                                    
#&gt; [225,] &quot;is_nonpositive_integer&quot;            &quot;is_nonzero_size&quot;                           
#&gt; [226,] &quot;is_not_nan&quot;                        &quot;is_ordered&quot;                                
#&gt; [227,] &quot;is_plain_type&quot;                     &quot;is_pos_definite&quot;                           
#&gt; [228,] &quot;is_positive&quot;                       &quot;is_rev_matrix&quot;                             
#&gt; [229,] &quot;is_scal_finite&quot;                    &quot;is_size_match&quot;                             
#&gt; [230,] &quot;is_square&quot;                         &quot;is_stan_scalar&quot;                            
#&gt; [231,] &quot;is_stan_scalar_or_eigen&quot;           &quot;is_string_convertible&quot;                     
#&gt; [232,] &quot;is_symmetric&quot;                      &quot;is_tuple&quot;                                  
#&gt; [233,] &quot;is_uninitialized&quot;                  &quot;is_unit_vector&quot;                            
#&gt; [234,] &quot;is_var&quot;                            &quot;is_var_and_matrix_types&quot;                   
#&gt; [235,] &quot;is_var_dense_dynamic&quot;              &quot;is_var_eigen&quot;                              
#&gt; [236,] &quot;is_var_matrix&quot;                     &quot;is_var_or_arithmetic&quot;                      
#&gt; [237,] &quot;is_vari&quot;                           &quot;is_vector&quot;                                 
#&gt; [238,] &quot;is_vector_like&quot;                    &quot;isfinite&quot;                                  
#&gt; [239,] &quot;isinf&quot;                             &quot;isnan&quot;                                     
#&gt; [240,] &quot;isnormal&quot;                          &quot;lambert_w&quot;                                 
#&gt; [241,] &quot;lb_constrain&quot;                      &quot;lb_free&quot;                                   
#&gt; [242,] &quot;lbeta&quot;                             &quot;ldexp&quot;                                     
#&gt; [243,] &quot;lgamma&quot;                            &quot;lgamma_stirling&quot;                           
#&gt; [244,] &quot;lgamma_stirling_diff&quot;              &quot;linspaced_array&quot;                           
#&gt; [245,] &quot;linspaced_int_array&quot;               &quot;linspaced_row_vector&quot;                      
#&gt; [246,] &quot;linspaced_vector&quot;                  &quot;lkj_corr_cholesky_log&quot;                     
#&gt; [247,] &quot;lkj_corr_cholesky_lpdf&quot;            &quot;lkj_corr_cholesky_rng&quot;                     
#&gt; [248,] &quot;lkj_corr_log&quot;                      &quot;lkj_corr_lpdf&quot;                             
#&gt; [249,] &quot;lkj_corr_rng&quot;                      &quot;lkj_cov_log&quot;                               
#&gt; [250,] &quot;lkj_cov_lpdf&quot;                      &quot;lmgamma&quot;                                   
#&gt; [251,] &quot;lmultiply&quot;                         &quot;log&quot;                                       
#&gt; [252,] &quot;log10&quot;                             &quot;log1m&quot;                                     
#&gt; [253,] &quot;log1m_exp&quot;                         &quot;log1m_inv_logit&quot;                           
#&gt; [254,] &quot;log1p&quot;                             &quot;log1p_exp&quot;                                 
#&gt; [255,] &quot;log2&quot;                              &quot;log_determinant&quot;                           
#&gt; [256,] &quot;log_determinant_ldlt&quot;              &quot;log_determinant_spd&quot;                       
#&gt; [257,] &quot;log_diff_exp&quot;                      &quot;log_falling_factorial&quot;                     
#&gt; [258,] &quot;log_inv_logit&quot;                     &quot;log_inv_logit_diff&quot;                        
#&gt; [259,] &quot;log_mix&quot;                           &quot;log_modified_bessel_first_kind&quot;            
#&gt; [260,] &quot;log_rising_factorial&quot;              &quot;log_softmax&quot;                               
#&gt; [261,] &quot;log_sum_exp&quot;                       &quot;log_sum_exp_signed&quot;                        
#&gt; [262,] &quot;logb&quot;                              &quot;logical_and&quot;                               
#&gt; [263,] &quot;logical_eq&quot;                        &quot;logical_gt&quot;                                
#&gt; [264,] &quot;logical_gte&quot;                       &quot;logical_lt&quot;                                
#&gt; [265,] &quot;logical_lte&quot;                       &quot;logical_negation&quot;                          
#&gt; [266,] &quot;logical_neq&quot;                       &quot;logical_or&quot;                                
#&gt; [267,] &quot;logistic_ccdf_log&quot;                 &quot;logistic_cdf&quot;                              
#&gt; [268,] &quot;logistic_cdf_log&quot;                  &quot;logistic_lccdf&quot;                            
#&gt; [269,] &quot;logistic_lcdf&quot;                     &quot;logistic_log&quot;                              
#&gt; [270,] &quot;logistic_lpdf&quot;                     &quot;logistic_rng&quot;                              
#&gt; [271,] &quot;logit&quot;                             &quot;loglogistic_cdf&quot;                           
#&gt; [272,] &quot;loglogistic_log&quot;                   &quot;loglogistic_lpdf&quot;                          
#&gt; [273,] &quot;loglogistic_rng&quot;                   &quot;lognormal_ccdf_log&quot;                        
#&gt; [274,] &quot;lognormal_cdf&quot;                     &quot;lognormal_cdf_log&quot;                         
#&gt; [275,] &quot;lognormal_lccdf&quot;                   &quot;lognormal_lcdf&quot;                            
#&gt; [276,] &quot;lognormal_log&quot;                     &quot;lognormal_lpdf&quot;                            
#&gt; [277,] &quot;lognormal_rng&quot;                     &quot;lub_constrain&quot;                             
#&gt; [278,] &quot;lub_free&quot;                          &quot;make_iter_name&quot;                            
#&gt; [279,] &quot;make_nu&quot;                           &quot;map_rect&quot;                                  
#&gt; [280,] &quot;map_rect_combine&quot;                  &quot;map_rect_concurrent&quot;                       
#&gt; [281,] &quot;map_rect_mpi&quot;                      &quot;map_rect_reduce&quot;                           
#&gt; [282,] &quot;matrix_exp&quot;                        &quot;matrix_exp_2x2&quot;                            
#&gt; [283,] &quot;matrix_exp_action_handler&quot;         &quot;matrix_exp_multiply&quot;                       
#&gt; [284,] &quot;matrix_exp_pade&quot;                   &quot;matrix_normal_prec_log&quot;                    
#&gt; [285,] &quot;matrix_normal_prec_lpdf&quot;           &quot;matrix_normal_prec_rng&quot;                    
#&gt; [286,] &quot;matrix_power&quot;                      &quot;max&quot;                                       
#&gt; [287,] &quot;max_size&quot;                          &quot;max_size_mvt&quot;                              
#&gt; [288,] &quot;mdivide_left&quot;                      &quot;mdivide_left_ldlt&quot;                         
#&gt; [289,] &quot;mdivide_left_spd&quot;                  &quot;mdivide_left_tri&quot;                          
#&gt; [290,] &quot;mdivide_left_tri_low&quot;              &quot;mdivide_right&quot;                             
#&gt; [291,] &quot;mdivide_right_ldlt&quot;                &quot;mdivide_right_spd&quot;                         
#&gt; [292,] &quot;mdivide_right_tri&quot;                 &quot;mdivide_right_tri_low&quot;                     
#&gt; [293,] &quot;mean&quot;                              &quot;min&quot;                                       
#&gt; [294,] &quot;minus&quot;                             &quot;modified_bessel_first_kind&quot;                
#&gt; [295,] &quot;modified_bessel_second_kind&quot;       &quot;modulus&quot;                                   
#&gt; [296,] &quot;mpi_cluster&quot;                       &quot;mpi_command&quot;                               
#&gt; [297,] &quot;mpi_distributed_apply&quot;             &quot;mpi_parallel_call&quot;                         
#&gt; [298,] &quot;multi_gp_cholesky_log&quot;             &quot;multi_gp_cholesky_lpdf&quot;                    
#&gt; [299,] &quot;multi_gp_log&quot;                      &quot;multi_gp_lpdf&quot;                             
#&gt; [300,] &quot;multi_normal_cholesky_log&quot;         &quot;multi_normal_cholesky_lpdf&quot;                
#&gt; [301,] &quot;multi_normal_cholesky_rng&quot;         &quot;multi_normal_log&quot;                          
#&gt; [302,] &quot;multi_normal_lpdf&quot;                 &quot;multi_normal_prec_log&quot;                     
#&gt; [303,] &quot;multi_normal_prec_lpdf&quot;            &quot;multi_normal_prec_rng&quot;                     
#&gt; [304,] &quot;multi_normal_rng&quot;                  &quot;multi_student_t_cholesky_lpdf&quot;             
#&gt; [305,] &quot;multi_student_t_cholesky_rng&quot;      &quot;multi_student_t_log&quot;                       
#&gt; [306,] &quot;multi_student_t_lpdf&quot;              &quot;multi_student_t_rng&quot;                       
#&gt; [307,] &quot;multinomial_log&quot;                   &quot;multinomial_logit_log&quot;                     
#&gt; [308,] &quot;multinomial_logit_lpmf&quot;            &quot;multinomial_logit_rng&quot;                     
#&gt; [309,] &quot;multinomial_lpmf&quot;                  &quot;multinomial_rng&quot;                           
#&gt; [310,] &quot;multiply&quot;                          &quot;multiply_log&quot;                              
#&gt; [311,] &quot;multiply_lower_tri_self_transpose&quot; &quot;neg_binomial_2_ccdf_log&quot;                   
#&gt; [312,] &quot;neg_binomial_2_cdf&quot;                &quot;neg_binomial_2_cdf_log&quot;                    
#&gt; [313,] &quot;neg_binomial_2_lccdf&quot;              &quot;neg_binomial_2_lcdf&quot;                       
#&gt; [314,] &quot;neg_binomial_2_log&quot;                &quot;neg_binomial_2_log_glm_log&quot;                
#&gt; [315,] &quot;neg_binomial_2_log_glm_lpmf&quot;       &quot;neg_binomial_2_log_log&quot;                    
#&gt; [316,] &quot;neg_binomial_2_log_lpmf&quot;           &quot;neg_binomial_2_log_rng&quot;                    
#&gt; [317,] &quot;neg_binomial_2_lpmf&quot;               &quot;neg_binomial_2_rng&quot;                        
#&gt; [318,] &quot;neg_binomial_ccdf_log&quot;             &quot;neg_binomial_cdf&quot;                          
#&gt; [319,] &quot;neg_binomial_cdf_log&quot;              &quot;neg_binomial_lccdf&quot;                        
#&gt; [320,] &quot;neg_binomial_lcdf&quot;                 &quot;neg_binomial_log&quot;                          
#&gt; [321,] &quot;neg_binomial_lpmf&quot;                 &quot;neg_binomial_rng&quot;                          
#&gt; [322,] &quot;norm&quot;                              &quot;norm1&quot;                                     
#&gt; [323,] &quot;norm2&quot;                             &quot;normal_ccdf_log&quot;                           
#&gt; [324,] &quot;normal_cdf&quot;                        &quot;normal_cdf_log&quot;                            
#&gt; [325,] &quot;normal_id_glm_log&quot;                 &quot;normal_id_glm_lpdf&quot;                        
#&gt; [326,] &quot;normal_lccdf&quot;                      &quot;normal_lcdf&quot;                               
#&gt; [327,] &quot;normal_log&quot;                        &quot;normal_lpdf&quot;                               
#&gt; [328,] &quot;normal_rng&quot;                        &quot;normal_sufficient_log&quot;                     
#&gt; [329,] &quot;normal_sufficient_lpdf&quot;            &quot;num_elements&quot;                              
#&gt; [330,] &quot;ode_ckrk&quot;                          &quot;ode_rk45&quot;                                  
#&gt; [331,] &quot;ode_store_sensitivities&quot;           &quot;offset_multiplier_constrain&quot;               
#&gt; [332,] &quot;offset_multiplier_free&quot;            &quot;one_hot_array&quot;                             
#&gt; [333,] &quot;one_hot_int_array&quot;                 &quot;one_hot_row_vector&quot;                        
#&gt; [334,] &quot;one_hot_vector&quot;                    &quot;ones_array&quot;                                
#&gt; [335,] &quot;ones_int_array&quot;                    &quot;ones_row_vector&quot;                           
#&gt; [336,] &quot;ones_vector&quot;                       &quot;operands_and_partials&quot;                     
#&gt; [337,] &quot;operator_addition&quot;                 &quot;operator_division&quot;                         
#&gt; [338,] &quot;operator_equal_equal&quot;              &quot;operator_minus&quot;                            
#&gt; [339,] &quot;operator_multiplication&quot;           &quot;operator_not_equal&quot;                        
#&gt; [340,] &quot;operator_plus&quot;                     &quot;operator_subtraction&quot;                      
#&gt; [341,] &quot;ordered_constrain&quot;                 &quot;ordered_free&quot;                              
#&gt; [342,] &quot;ordered_logistic_glm_lpmf&quot;         &quot;ordered_logistic_log&quot;                      
#&gt; [343,] &quot;ordered_logistic_lpmf&quot;             &quot;ordered_logistic_rng&quot;                      
#&gt; [344,] &quot;ordered_probit_log&quot;                &quot;ordered_probit_lpmf&quot;                       
#&gt; [345,] &quot;ordered_probit_rng&quot;                &quot;out_of_range&quot;                              
#&gt; [346,] &quot;owens_t&quot;                           &quot;pareto_ccdf_log&quot;                           
#&gt; [347,] &quot;pareto_cdf&quot;                        &quot;pareto_cdf_log&quot;                            
#&gt; [348,] &quot;pareto_lccdf&quot;                      &quot;pareto_lcdf&quot;                               
#&gt; [349,] &quot;pareto_log&quot;                        &quot;pareto_lpdf&quot;                               
#&gt; [350,] &quot;pareto_rng&quot;                        &quot;pareto_type_2_ccdf_log&quot;                    
#&gt; [351,] &quot;pareto_type_2_cdf&quot;                 &quot;pareto_type_2_cdf_log&quot;                     
#&gt; [352,] &quot;pareto_type_2_lccdf&quot;               &quot;pareto_type_2_lcdf&quot;                        
#&gt; [353,] &quot;pareto_type_2_log&quot;                 &quot;pareto_type_2_lpdf&quot;                        
#&gt; [354,] &quot;pareto_type_2_rng&quot;                 &quot;partials_propagator&quot;                       
#&gt; [355,] &quot;partials_return_type&quot;              &quot;partials_type&quot;                             
#&gt; [356,] &quot;plain_type&quot;                        &quot;plus&quot;                                      
#&gt; [357,] &quot;poisson_binomial_ccdf_log&quot;         &quot;poisson_binomial_cdf&quot;                      
#&gt; [358,] &quot;poisson_binomial_cdf_log&quot;          &quot;poisson_binomial_lccdf&quot;                    
#&gt; [359,] &quot;poisson_binomial_lcdf&quot;             &quot;poisson_binomial_log&quot;                      
#&gt; [360,] &quot;poisson_binomial_log_probs&quot;        &quot;poisson_binomial_lpmf&quot;                     
#&gt; [361,] &quot;poisson_binomial_rng&quot;              &quot;poisson_ccdf_log&quot;                          
#&gt; [362,] &quot;poisson_cdf&quot;                       &quot;poisson_cdf_log&quot;                           
#&gt; [363,] &quot;poisson_lccdf&quot;                     &quot;poisson_lcdf&quot;                              
#&gt; [364,] &quot;poisson_log&quot;                       &quot;poisson_log_glm_log&quot;                       
#&gt; [365,] &quot;poisson_log_glm_lpmf&quot;              &quot;poisson_log_log&quot;                           
#&gt; [366,] &quot;poisson_log_lpmf&quot;                  &quot;poisson_log_rng&quot;                           
#&gt; [367,] &quot;poisson_lpmf&quot;                      &quot;poisson_rng&quot;                               
#&gt; [368,] &quot;polar&quot;                             &quot;positive_constrain&quot;                        
#&gt; [369,] &quot;positive_free&quot;                     &quot;positive_ordered_constrain&quot;                
#&gt; [370,] &quot;positive_ordered_free&quot;             &quot;possibly_sum&quot;                              
#&gt; [371,] &quot;pow&quot;                               &quot;primitive_value&quot;                           
#&gt; [372,] &quot;prob_constrain&quot;                    &quot;prob_free&quot;                                 
#&gt; [373,] &quot;prod&quot;                              &quot;proj&quot;                                      
#&gt; [374,] &quot;promote_args&quot;                      &quot;promote_elements&quot;                          
#&gt; [375,] &quot;promote_scalar&quot;                    &quot;promote_scalar_type&quot;                       
#&gt; [376,] &quot;pseudo_eigenvalues&quot;                &quot;pseudo_eigenvectors&quot;                       
#&gt; [377,] &quot;qr&quot;                                &quot;qr_Q&quot;                                      
#&gt; [378,] &quot;qr_R&quot;                              &quot;qr_thin&quot;                                   
#&gt; [379,] &quot;qr_thin_Q&quot;                         &quot;qr_thin_R&quot;                                 
#&gt; [380,] &quot;quad_form&quot;                         &quot;quad_form_diag&quot;                            
#&gt; [381,] &quot;quad_form_sym&quot;                     &quot;quantile&quot;                                  
#&gt; [382,] &quot;rank&quot;                              &quot;rayleigh_ccdf_log&quot;                         
#&gt; [383,] &quot;rayleigh_cdf&quot;                      &quot;rayleigh_cdf_log&quot;                          
#&gt; [384,] &quot;rayleigh_lccdf&quot;                    &quot;rayleigh_lcdf&quot;                             
#&gt; [385,] &quot;rayleigh_log&quot;                      &quot;rayleigh_lpdf&quot;                             
#&gt; [386,] &quot;rayleigh_rng&quot;                      &quot;read_corr_L&quot;                               
#&gt; [387,] &quot;read_corr_matrix&quot;                  &quot;read_cov_L&quot;                                
#&gt; [388,] &quot;read_cov_matrix&quot;                   &quot;real&quot;                                      
#&gt; [389,] &quot;reduce_sum&quot;                        &quot;reduce_sum_static&quot;                         
#&gt; [390,] &quot;ref_type&quot;                          &quot;rep_array&quot;                                 
#&gt; [391,] &quot;rep_matrix&quot;                        &quot;rep_row_vector&quot;                            
#&gt; [392,] &quot;rep_vector&quot;                        &quot;require_generics&quot;                          
#&gt; [393,] &quot;require_helpers&quot;                   &quot;resize&quot;                                    
#&gt; [394,] &quot;return_type&quot;                       &quot;reverse&quot;                                   
#&gt; [395,] &quot;rising_factorial&quot;                  &quot;round&quot;                                     
#&gt; [396,] &quot;row&quot;                               &quot;rows&quot;                                      
#&gt; [397,] &quot;rows_dot_product&quot;                  &quot;rows_dot_self&quot;                             
#&gt; [398,] &quot;scalar_seq_view&quot;                   &quot;scalar_type&quot;                               
#&gt; [399,] &quot;scalar_type_pre&quot;                   &quot;scalbn&quot;                                    
#&gt; [400,] &quot;scale_matrix_exp_multiply&quot;         &quot;scaled_add&quot;                                
#&gt; [401,] &quot;scaled_inv_chi_square_ccdf_log&quot;    &quot;scaled_inv_chi_square_cdf&quot;                 
#&gt; [402,] &quot;scaled_inv_chi_square_cdf_log&quot;     &quot;scaled_inv_chi_square_lccdf&quot;               
#&gt; [403,] &quot;scaled_inv_chi_square_lcdf&quot;        &quot;scaled_inv_chi_square_log&quot;                 
#&gt; [404,] &quot;scaled_inv_chi_square_lpdf&quot;        &quot;scaled_inv_chi_square_rng&quot;                 
#&gt; [405,] &quot;sd&quot;                                &quot;segment&quot;                                   
#&gt; [406,] &quot;select&quot;                            &quot;seq_view&quot;                                  
#&gt; [407,] &quot;sign&quot;                              &quot;signbit&quot;                                   
#&gt; [408,] &quot;simplex_constrain&quot;                 &quot;simplex_free&quot;                              
#&gt; [409,] &quot;sin&quot;                               &quot;singular_values&quot;                           
#&gt; [410,] &quot;sinh&quot;                              &quot;size&quot;                                      
#&gt; [411,] &quot;size_mvt&quot;                          &quot;size_zero&quot;                                 
#&gt; [412,] &quot;skew_double_exponential_ccdf_log&quot;  &quot;skew_double_exponential_cdf&quot;               
#&gt; [413,] &quot;skew_double_exponential_cdf_log&quot;   &quot;skew_double_exponential_lccdf&quot;             
#&gt; [414,] &quot;skew_double_exponential_lcdf&quot;      &quot;skew_double_exponential_log&quot;               
#&gt; [415,] &quot;skew_double_exponential_lpdf&quot;      &quot;skew_double_exponential_rng&quot;               
#&gt; [416,] &quot;skew_normal_ccdf_log&quot;              &quot;skew_normal_cdf&quot;                           
#&gt; [417,] &quot;skew_normal_cdf_log&quot;               &quot;skew_normal_lccdf&quot;                         
#&gt; [418,] &quot;skew_normal_lcdf&quot;                  &quot;skew_normal_log&quot;                           
#&gt; [419,] &quot;skew_normal_lpdf&quot;                  &quot;skew_normal_rng&quot;                           
#&gt; [420,] &quot;softmax&quot;                           &quot;sort_asc&quot;                                  
#&gt; [421,] &quot;sort_desc&quot;                         &quot;sort_indices&quot;                              
#&gt; [422,] &quot;sort_indices_asc&quot;                  &quot;sort_indices_desc&quot;                         
#&gt; [423,] &quot;sqrt&quot;                              &quot;square&quot;                                    
#&gt; [424,] &quot;squared_distance&quot;                  &quot;stan_print&quot;                                
#&gt; [425,] &quot;static_select&quot;                     &quot;std_normal_ccdf_log&quot;                       
#&gt; [426,] &quot;std_normal_cdf&quot;                    &quot;std_normal_cdf_log&quot;                        
#&gt; [427,] &quot;std_normal_lccdf&quot;                  &quot;std_normal_lcdf&quot;                           
#&gt; [428,] &quot;std_normal_log&quot;                    &quot;std_normal_log_qf&quot;                         
#&gt; [429,] &quot;std_normal_lpdf&quot;                   &quot;std_normal_rng&quot;                            
#&gt; [430,] &quot;step&quot;                              &quot;student_t_ccdf_log&quot;                        
#&gt; [431,] &quot;student_t_cdf&quot;                     &quot;student_t_cdf_log&quot;                         
#&gt; [432,] &quot;student_t_lccdf&quot;                   &quot;student_t_lcdf&quot;                            
#&gt; [433,] &quot;student_t_log&quot;                     &quot;student_t_lpdf&quot;                            
#&gt; [434,] &quot;student_t_rng&quot;                     &quot;sub_col&quot;                                   
#&gt; [435,] &quot;sub_row&quot;                           &quot;subtract&quot;                                  
#&gt; [436,] &quot;sum&quot;                               &quot;svd&quot;                                       
#&gt; [437,] &quot;svd_U&quot;                             &quot;svd_V&quot;                                     
#&gt; [438,] &quot;symmetrize_from_lower_tri&quot;         &quot;symmetrize_from_upper_tri&quot;                 
#&gt; [439,] &quot;system_error&quot;                      &quot;tail&quot;                                      
#&gt; [440,] &quot;tan&quot;                               &quot;tanh&quot;                                      
#&gt; [441,] &quot;tcrossprod&quot;                        &quot;tgamma&quot;                                    
#&gt; [442,] &quot;throw_domain_error&quot;                &quot;throw_domain_error_mat&quot;                    
#&gt; [443,] &quot;throw_domain_error_vec&quot;            &quot;to_array_1d&quot;                               
#&gt; [444,] &quot;to_array_2d&quot;                       &quot;to_complex&quot;                                
#&gt; [445,] &quot;to_int&quot;                            &quot;to_matrix&quot;                                 
#&gt; [446,] &quot;to_ref&quot;                            &quot;to_row_vector&quot;                             
#&gt; [447,] &quot;to_vector&quot;                         &quot;trace&quot;                                     
#&gt; [448,] &quot;trace_gen_inv_quad_form_ldlt&quot;      &quot;trace_gen_quad_form&quot;                       
#&gt; [449,] &quot;trace_inv_quad_form_ldlt&quot;          &quot;trace_quad_form&quot;                           
#&gt; [450,] &quot;transpose&quot;                         &quot;trigamma&quot;                                  
#&gt; [451,] &quot;trunc&quot;                             &quot;typedefs&quot;                                  
#&gt; [452,] &quot;ub_constrain&quot;                      &quot;ub_free&quot;                                   
#&gt; [453,] &quot;uniform_ccdf_log&quot;                  &quot;uniform_cdf&quot;                               
#&gt; [454,] &quot;uniform_cdf_log&quot;                   &quot;uniform_lccdf&quot;                             
#&gt; [455,] &quot;uniform_lcdf&quot;                      &quot;uniform_log&quot;                               
#&gt; [456,] &quot;uniform_lpdf&quot;                      &quot;uniform_rng&quot;                               
#&gt; [457,] &quot;uniform_simplex&quot;                   &quot;unit_vector_constrain&quot;                     
#&gt; [458,] &quot;unit_vector_free&quot;                  &quot;unitspaced_array&quot;                          
#&gt; [459,] &quot;validate_non_negative_index&quot;       &quot;validate_positive_index&quot;                   
#&gt; [460,] &quot;validate_unit_vector_index&quot;        &quot;value_of&quot;                                  
#&gt; [461,] &quot;value_of_rec&quot;                      &quot;value_type&quot;                                
#&gt; [462,] &quot;variance&quot;                          &quot;vec_concat&quot;                                
#&gt; [463,] &quot;vector_seq_view&quot;                   &quot;void_t&quot;                                    
#&gt; [464,] &quot;von_mises_ccdf_log&quot;                &quot;von_mises_cdf&quot;                             
#&gt; [465,] &quot;von_mises_cdf_log&quot;                 &quot;von_mises_lccdf&quot;                           
#&gt; [466,] &quot;von_mises_lcdf&quot;                    &quot;von_mises_log&quot;                             
#&gt; [467,] &quot;von_mises_lpdf&quot;                    &quot;von_mises_rng&quot;                             
#&gt; [468,] &quot;weibull_ccdf_log&quot;                  &quot;weibull_cdf&quot;                               
#&gt; [469,] &quot;weibull_cdf_log&quot;                   &quot;weibull_lccdf&quot;                             
#&gt; [470,] &quot;weibull_lcdf&quot;                      &quot;weibull_log&quot;                               
#&gt; [471,] &quot;weibull_lpdf&quot;                      &quot;weibull_rng&quot;                               
#&gt; [472,] &quot;welford_covar_estimator&quot;           &quot;welford_var_estimator&quot;                     
#&gt; [473,] &quot;wiener_log&quot;                        &quot;wiener_lpdf&quot;                               
#&gt; [474,] &quot;wishart_cholesky_lpdf&quot;             &quot;wishart_cholesky_rng&quot;                      
#&gt; [475,] &quot;wishart_log&quot;                       &quot;wishart_lpdf&quot;                              
#&gt; [476,] &quot;wishart_rng&quot;                       &quot;zeros_array&quot;                               
#&gt; [477,] &quot;zeros_int_array&quot;                   &quot;zeros_row_vector&quot;                          
#&gt; [478,] &quot;zeros_vector&quot;                      &quot;&quot;</code></pre>
</div>
<div id="using-higher-order-functions-in-the-stanheaders-package" class="section level1">
<h1>Using Higher-Order Functions in the <strong>StanHeaders</strong> Package</h1>
<p>This section will demonstrate how to use some of the C++ functions in the <strong>StanHeaders</strong> package whose first argument is another C++ function, in which case the <code>stanFunction</code> in the previous section will not work and you have to write your own C++.</p>
<div id="derivatives-and-minimization" class="section level2">
<h2>Derivatives and Minimization</h2>
<p>The following is a toy example of using the Stan Math library via <code>Rcpp::sourceCpp</code>: to minimize the function <span class="math display">\[\left(\mathbf{x} - \mathbf{a}\right)^\top \left(\mathbf{x} - \mathbf{a}\right)\]</span> which has a global minimum when <span class="math inline">\(\mathbf{x} = \mathbf{a}\)</span>. To find this minimum with autodifferentiation, we need to define the objective function. Then, its gradient with respect to <span class="math inline">\(\mathbf{x}\)</span>, which we know is <span class="math inline">\(2\left(\mathbf{x} - \mathbf{a}\right)\)</span> in this case, can be calculated by autodifferentiation. At the optimum (or on the way to the optimum), we might want to evaluate the Hessian matrix, which we know is <span class="math inline">\(2\mathbf{I}\)</span>, but would need an additional function to evaluate it via autodifferentiation. Finally, one could reconceptualize the problem as solving a homogeneous system of equations where the gradient is set equal to a vector of zeros. The <code>stan::math::algebra_solver</code> function can solve such a system using autodifferentiation to obtain the Jacobian, which we know to be the identity matrix in this case.</p>
<div class="sourceCode" id="cb6"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb6-1" title="1"><span class="kw">Sys.setenv</span>(<span class="dt">PKG_CXXFLAGS =</span> StanHeaders<span class="op">:::</span><span class="kw">CxxFlags</span>(<span class="dt">as_character =</span> <span class="ot">TRUE</span>))</a>
<a class="sourceLine" id="cb6-2" title="2">SH &lt;-<span class="st"> </span><span class="kw">system.file</span>(<span class="kw">ifelse</span>(.Platform<span class="op">$</span>OS.type <span class="op">==</span><span class="st"> &quot;windows&quot;</span>, <span class="st">&quot;libs&quot;</span>, <span class="st">&quot;lib&quot;</span>), .Platform<span class="op">$</span>r_arch,</a>
<a class="sourceLine" id="cb6-3" title="3">                  <span class="dt">package =</span> <span class="st">&quot;StanHeaders&quot;</span>, <span class="dt">mustWork =</span> <span class="ot">TRUE</span>)</a>
<a class="sourceLine" id="cb6-4" title="4"><span class="kw">Sys.setenv</span>(<span class="dt">PKG_LIBS =</span> <span class="kw">paste0</span>(StanHeaders<span class="op">:::</span><span class="kw">LdFlags</span>(<span class="dt">as_character =</span> <span class="ot">TRUE</span>),</a>
<a class="sourceLine" id="cb6-5" title="5">                             <span class="st">&quot; -L&quot;</span>, <span class="kw">shQuote</span>(SH), <span class="st">&quot; -lStanHeaders&quot;</span>))</a>
<a class="sourceLine" id="cb6-6" title="6">                            </a></code></pre></div>
<p>Here is C++ code that does all of the above, except for the part of finding the optimum, which is done using the R function <code>optim</code> below.</p>
<div class="sourceCode" id="cb7"><pre class="sourceCode cpp"><code class="sourceCode cpp"><a class="sourceLine" id="cb7-1" title="1"><span class="co">// [[Rcpp::depends(BH)]]</span></a>
<a class="sourceLine" id="cb7-2" title="2"><span class="co">// [[Rcpp::depends(RcppEigen)]]</span></a>
<a class="sourceLine" id="cb7-3" title="3"><span class="co">// [[Rcpp::depends(RcppParallel)]]</span></a>
<a class="sourceLine" id="cb7-4" title="4"><span class="co">// [[Rcpp::depends(StanHeaders)]]</span></a>
<a class="sourceLine" id="cb7-5" title="5"><span class="pp">#include </span><span class="im">&lt;stan/math/mix.hpp&gt;</span><span class="pp"> </span><span class="co">// stuff from mix/ must come first</span></a>
<a class="sourceLine" id="cb7-6" title="6"><span class="pp">#include </span><span class="im">&lt;stan/math.hpp&gt;</span><span class="pp">     </span><span class="co">// finally pull in everything from rev/ and prim/</span></a>
<a class="sourceLine" id="cb7-7" title="7"><span class="pp">#include </span><span class="im">&lt;Rcpp.h&gt;</span></a>
<a class="sourceLine" id="cb7-8" title="8"><span class="pp">#include </span><span class="im">&lt;RcppEigen.h&gt;</span><span class="pp">       </span><span class="co">// do this AFTER including stuff from stan/math</span></a>
<a class="sourceLine" id="cb7-9" title="9"></a>
<a class="sourceLine" id="cb7-10" title="10"><span class="co">// [[Rcpp::plugins(cpp17)]]</span></a>
<a class="sourceLine" id="cb7-11" title="11"></a>
<a class="sourceLine" id="cb7-12" title="12"><span class="co">/* Objective function */</span></a>
<a class="sourceLine" id="cb7-13" title="13"></a>
<a class="sourceLine" id="cb7-14" title="14"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb7-15" title="15"><span class="kw">auto</span> f(Eigen::VectorXd x, Eigen::VectorXd a) { <span class="co">// objective function in doubles</span></a>
<a class="sourceLine" id="cb7-16" title="16">  <span class="kw">using</span> stan::math::dot_self;                  <span class="co">// dot_self() is a dot product with self</span></a>
<a class="sourceLine" id="cb7-17" title="17">  <span class="cf">return</span> dot_self( (x - a) );</a>
<a class="sourceLine" id="cb7-18" title="18">}</a>
<a class="sourceLine" id="cb7-19" title="19"></a>
<a class="sourceLine" id="cb7-20" title="20"><span class="co">/* Gradient */</span></a>
<a class="sourceLine" id="cb7-21" title="21"></a>
<a class="sourceLine" id="cb7-22" title="22"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb7-23" title="23"><span class="kw">auto</span> g(Eigen::VectorXd x, Eigen::VectorXd a) {  <span class="co">// gradient by AD using Stan</span></a>
<a class="sourceLine" id="cb7-24" title="24">  <span class="dt">double</span> fx;</a>
<a class="sourceLine" id="cb7-25" title="25">  Eigen::VectorXd grad_fx;</a>
<a class="sourceLine" id="cb7-26" title="26">  <span class="kw">using</span> stan::math::dot_self;</a>
<a class="sourceLine" id="cb7-27" title="27">  stan::math::gradient([&amp;a](<span class="kw">auto</span> x) { <span class="cf">return</span> dot_self( (x - a) ); },</a>
<a class="sourceLine" id="cb7-28" title="28">                       x, fx, grad_fx);</a>
<a class="sourceLine" id="cb7-29" title="29">  <span class="cf">return</span> grad_fx;</a>
<a class="sourceLine" id="cb7-30" title="30">}</a>
<a class="sourceLine" id="cb7-31" title="31"></a>
<a class="sourceLine" id="cb7-32" title="32"><span class="co">/* Hessian */</span></a>
<a class="sourceLine" id="cb7-33" title="33"></a>
<a class="sourceLine" id="cb7-34" title="34"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb7-35" title="35"><span class="kw">auto</span> H(Eigen::VectorXd x, Eigen::VectorXd a) { <span class="co">// Hessian by AD using Stan</span></a>
<a class="sourceLine" id="cb7-36" title="36">  <span class="dt">double</span> fx;</a>
<a class="sourceLine" id="cb7-37" title="37">  Eigen::VectorXd grad_fx;</a>
<a class="sourceLine" id="cb7-38" title="38">  Eigen::MatrixXd H;</a>
<a class="sourceLine" id="cb7-39" title="39">  <span class="kw">using</span> stan::math::dot_self;</a>
<a class="sourceLine" id="cb7-40" title="40">  stan::math::hessian([&amp;a](<span class="kw">auto</span> x) { <span class="cf">return</span> dot_self(x - a); },</a>
<a class="sourceLine" id="cb7-41" title="41">                      x, fx, grad_fx, H);</a>
<a class="sourceLine" id="cb7-42" title="42">  <span class="cf">return</span> H;</a>
<a class="sourceLine" id="cb7-43" title="43">}</a>
<a class="sourceLine" id="cb7-44" title="44"></a>
<a class="sourceLine" id="cb7-45" title="45"><span class="co">/* Jacobian */</span></a>
<a class="sourceLine" id="cb7-46" title="46"></a>
<a class="sourceLine" id="cb7-47" title="47"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb7-48" title="48"><span class="kw">auto</span> J(Eigen::VectorXd x, Eigen::VectorXd a) { <span class="co">// not actually used</span></a>
<a class="sourceLine" id="cb7-49" title="49">  Eigen::VectorXd fx;</a>
<a class="sourceLine" id="cb7-50" title="50">  Eigen::MatrixXd J;</a>
<a class="sourceLine" id="cb7-51" title="51">  <span class="kw">using</span> stan::math::dot_self;</a>
<a class="sourceLine" id="cb7-52" title="52">  stan::math::jacobian([&amp;a](<span class="kw">auto</span> x) {</a>
<a class="sourceLine" id="cb7-53" title="53">    <span class="cf">return</span> (<span class="dv">2</span> * (x - a));</a>
<a class="sourceLine" id="cb7-54" title="54">  }, x, fx, J);</a>
<a class="sourceLine" id="cb7-55" title="55">  <span class="cf">return</span> J;</a>
<a class="sourceLine" id="cb7-56" title="56">}</a>
<a class="sourceLine" id="cb7-57" title="57"></a>
<a class="sourceLine" id="cb7-58" title="58"><span class="kw">struct</span> equations_functor {</a>
<a class="sourceLine" id="cb7-59" title="59">  <span class="kw">template</span> &lt;<span class="kw">typename</span> T0, <span class="kw">typename</span> T1&gt;</a>
<a class="sourceLine" id="cb7-60" title="60">  <span class="kw">inline</span> Eigen::Matrix&lt;T0, Eigen::Dynamic, <span class="dv">1</span>&gt;</a>
<a class="sourceLine" id="cb7-61" title="61">  <span class="kw">operator</span>()(<span class="at">const</span> Eigen::Matrix&lt;T0, Eigen::Dynamic, <span class="dv">1</span>&gt;&amp; x,</a>
<a class="sourceLine" id="cb7-62" title="62">             <span class="at">const</span> Eigen::Matrix&lt;T1, Eigen::Dynamic, <span class="dv">1</span>&gt;&amp; theta,</a>
<a class="sourceLine" id="cb7-63" title="63">             <span class="at">const</span> <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt;&amp; x_r, <span class="at">const</span> <span class="bu">std::</span>vector&lt;<span class="dt">int</span>&gt;&amp; x_i,</a>
<a class="sourceLine" id="cb7-64" title="64">             <span class="bu">std::</span>ostream* pstream__) <span class="at">const</span> {</a>
<a class="sourceLine" id="cb7-65" title="65">    <span class="cf">return</span> <span class="dv">2</span> * (x - stan::math::to_vector(x_r));</a>
<a class="sourceLine" id="cb7-66" title="66">  }</a>
<a class="sourceLine" id="cb7-67" title="67">};</a>
<a class="sourceLine" id="cb7-68" title="68"></a>
<a class="sourceLine" id="cb7-69" title="69"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb7-70" title="70"><span class="kw">auto</span> solution(Eigen::VectorXd a, Eigen::VectorXd guess) {</a>
<a class="sourceLine" id="cb7-71" title="71">  Eigen::VectorXd theta;</a>
<a class="sourceLine" id="cb7-72" title="72">  <span class="kw">auto</span> x_r = stan::math::to_array_1d(a);</a>
<a class="sourceLine" id="cb7-73" title="73">  equations_functor f;</a>
<a class="sourceLine" id="cb7-74" title="74">  <span class="kw">auto</span> x = stan::math::algebra_solver(f, guess, theta, x_r, {});</a>
<a class="sourceLine" id="cb7-75" title="75">  <span class="cf">return</span> x;</a>
<a class="sourceLine" id="cb7-76" title="76">}</a></code></pre></div>
<p>In this compiled RMarkdown document, the <strong>knitr</strong> package has exported functions <code>f</code>, <code>g</code>, <code>H</code>, <code>J</code> and <code>solution</code> (but not <code>equations_functor</code>) to R’s global environment using the <code>sourceCpp</code> function in the <strong>Rcpp</strong> package, so that they can now be called from R. Here we find the optimum starting from a random point in three dimensions:</p>
<div class="sourceCode" id="cb8"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb8-1" title="1">x &lt;-<span class="st"> </span><span class="kw">optim</span>(<span class="kw">rnorm</span>(<span class="dv">3</span>), <span class="dt">fn =</span> f, <span class="dt">gr =</span> g, <span class="dt">a =</span> <span class="dv">1</span><span class="op">:</span><span class="dv">3</span>, <span class="dt">method =</span> <span class="st">&quot;BFGS&quot;</span>, <span class="dt">hessian =</span> <span class="ot">TRUE</span>)</a>
<a class="sourceLine" id="cb8-2" title="2">x<span class="op">$</span>par</a>
<a class="sourceLine" id="cb8-3" title="3"><span class="co">#&gt; [1] 1 2 3</span></a>
<a class="sourceLine" id="cb8-4" title="4">x<span class="op">$</span>hessian</a>
<a class="sourceLine" id="cb8-5" title="5"><span class="co">#&gt;      [,1] [,2] [,3]</span></a>
<a class="sourceLine" id="cb8-6" title="6"><span class="co">#&gt; [1,]    2    0    0</span></a>
<a class="sourceLine" id="cb8-7" title="7"><span class="co">#&gt; [2,]    0    2    0</span></a>
<a class="sourceLine" id="cb8-8" title="8"><span class="co">#&gt; [3,]    0    0    2</span></a>
<a class="sourceLine" id="cb8-9" title="9"><span class="kw">H</span>(x<span class="op">$</span>par, <span class="dt">a =</span> <span class="dv">1</span><span class="op">:</span><span class="dv">3</span>)</a>
<a class="sourceLine" id="cb8-10" title="10"><span class="co">#&gt;      [,1] [,2] [,3]</span></a>
<a class="sourceLine" id="cb8-11" title="11"><span class="co">#&gt; [1,]    2    0    0</span></a>
<a class="sourceLine" id="cb8-12" title="12"><span class="co">#&gt; [2,]    0    2    0</span></a>
<a class="sourceLine" id="cb8-13" title="13"><span class="co">#&gt; [3,]    0    0    2</span></a>
<a class="sourceLine" id="cb8-14" title="14"><span class="kw">J</span>(x<span class="op">$</span>par, <span class="dt">a =</span> <span class="dv">1</span><span class="op">:</span><span class="dv">3</span>)</a>
<a class="sourceLine" id="cb8-15" title="15"><span class="co">#&gt;      [,1] [,2] [,3]</span></a>
<a class="sourceLine" id="cb8-16" title="16"><span class="co">#&gt; [1,]    2    0    0</span></a>
<a class="sourceLine" id="cb8-17" title="17"><span class="co">#&gt; [2,]    0    2    0</span></a>
<a class="sourceLine" id="cb8-18" title="18"><span class="co">#&gt; [3,]    0    0    2</span></a>
<a class="sourceLine" id="cb8-19" title="19"><span class="kw">solution</span>(<span class="dt">a =</span> <span class="dv">1</span><span class="op">:</span><span class="dv">3</span>, <span class="dt">guess =</span> <span class="kw">rnorm</span>(<span class="dv">3</span>))</a>
<a class="sourceLine" id="cb8-20" title="20"><span class="co">#&gt; [1] 1 2 3</span></a></code></pre></div>
</div>
</div>
<div id="integrals-and-ordinary-differential-equations" class="section level1">
<h1>Integrals and Ordinary Differential Equations</h1>
<p>The Stan Math library can do one-dimensional numerical integration and can solve stiff and non-stiff systems of differential equations, such as the harmonic oscillator example below. Solving stiff systems utilizes the CVODES library, which is included in <strong>StanHeaders</strong>.</p>
<div class="sourceCode" id="cb9"><pre class="sourceCode cpp"><code class="sourceCode cpp"><a class="sourceLine" id="cb9-1" title="1"><span class="co">// [[Rcpp::depends(BH)]]</span></a>
<a class="sourceLine" id="cb9-2" title="2"><span class="co">// [[Rcpp::depends(RcppEigen)]]</span></a>
<a class="sourceLine" id="cb9-3" title="3"><span class="co">// [[Rcpp::depends(RcppParallel)]]</span></a>
<a class="sourceLine" id="cb9-4" title="4"><span class="co">// [[Rcpp::depends(StanHeaders)]]</span></a>
<a class="sourceLine" id="cb9-5" title="5"><span class="pp">#include </span><span class="im">&lt;stan/math.hpp&gt;</span><span class="pp">                         </span><span class="co">// pulls in everything from rev/ and prim/</span></a>
<a class="sourceLine" id="cb9-6" title="6"><span class="pp">#include </span><span class="im">&lt;Rcpp.h&gt;</span></a>
<a class="sourceLine" id="cb9-7" title="7"><span class="pp">#include </span><span class="im">&lt;RcppEigen.h&gt;</span><span class="pp">                           </span><span class="co">// do this AFTER including stan/math</span></a>
<a class="sourceLine" id="cb9-8" title="8"></a>
<a class="sourceLine" id="cb9-9" title="9"><span class="co">// [[Rcpp::plugins(cpp17)]]</span></a>
<a class="sourceLine" id="cb9-10" title="10"></a>
<a class="sourceLine" id="cb9-11" title="11"><span class="co">/* Definite integrals */</span></a>
<a class="sourceLine" id="cb9-12" title="12"></a>
<a class="sourceLine" id="cb9-13" title="13"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb9-14" title="14"><span class="dt">double</span> Cauchy(<span class="dt">double</span> scale) {</a>
<a class="sourceLine" id="cb9-15" title="15">  <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt; theta;</a>
<a class="sourceLine" id="cb9-16" title="16">  <span class="kw">auto</span> half = stan::math::integrate_1d([](<span class="kw">auto</span> x, <span class="kw">auto</span> xc, <span class="kw">auto</span> theta,</a>
<a class="sourceLine" id="cb9-17" title="17">                                          <span class="kw">auto</span> x_r, <span class="kw">auto</span> x_i, <span class="kw">auto</span> msgs) {</a>
<a class="sourceLine" id="cb9-18" title="18">    <span class="cf">return</span> exp(stan::math::cauchy_lpdf(x, <span class="dv">0</span>, x_r[<span class="dv">0</span>]));</a>
<a class="sourceLine" id="cb9-19" title="19">  }, -scale, scale, theta, {scale}, {}, <span class="kw">nullptr</span>, <span class="fl">1e-7</span>);</a>
<a class="sourceLine" id="cb9-20" title="20">  <span class="cf">return</span> half * <span class="dv">2</span>; <span class="co">// should equal 1 for any positive scale</span></a>
<a class="sourceLine" id="cb9-21" title="21">}</a>
<a class="sourceLine" id="cb9-22" title="22"></a>
<a class="sourceLine" id="cb9-23" title="23"><span class="co">/* Ordinary Differential Equations */</span></a>
<a class="sourceLine" id="cb9-24" title="24"></a>
<a class="sourceLine" id="cb9-25" title="25"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb9-26" title="26"><span class="kw">auto</span> nonstiff(Eigen::MatrixXd A, Eigen::VectorXd y0) {</a>
<a class="sourceLine" id="cb9-27" title="27">  <span class="kw">using</span> stan::math::integrate_ode_rk45;</a>
<a class="sourceLine" id="cb9-28" title="28">  <span class="kw">using</span> stan::math::to_vector;</a>
<a class="sourceLine" id="cb9-29" title="29">  <span class="kw">using</span> stan::math::to_array_1d;</a>
<a class="sourceLine" id="cb9-30" title="30">  <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt; theta;</a>
<a class="sourceLine" id="cb9-31" title="31">  <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt; times = {<span class="dv">1</span>, <span class="dv">2</span>};</a>
<a class="sourceLine" id="cb9-32" title="32">  <span class="kw">auto</span> y = integrate_ode_rk45([&amp;A](<span class="kw">auto</span> t, <span class="kw">auto</span> y, </a>
<a class="sourceLine" id="cb9-33" title="33">                                   <span class="kw">auto</span> theta, <span class="kw">auto</span> x_r, <span class="kw">auto</span> x_i, <span class="bu">std::</span>ostream *msgs) {</a>
<a class="sourceLine" id="cb9-34" title="34">    <span class="cf">return</span> to_array_1d( (A * to_vector(y)).eval() );</a>
<a class="sourceLine" id="cb9-35" title="35">  }, to_array_1d(y0), <span class="dv">0</span>, times, theta, {}, {});</a>
<a class="sourceLine" id="cb9-36" title="36">  Eigen::VectorXd truth = stan::math::matrix_exp(A) * y0;</a>
<a class="sourceLine" id="cb9-37" title="37">  <span class="cf">return</span> (to_vector(y[<span class="dv">0</span>]) - truth).eval(); <span class="co">// should be &quot;zero&quot;</span></a>
<a class="sourceLine" id="cb9-38" title="38">}</a>
<a class="sourceLine" id="cb9-39" title="39"></a>
<a class="sourceLine" id="cb9-40" title="40"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb9-41" title="41"><span class="kw">auto</span> stiff(Eigen::MatrixXd A, Eigen::VectorXd y0) { <span class="co">// not actually stiff</span></a>
<a class="sourceLine" id="cb9-42" title="42">  <span class="kw">using</span> stan::math::integrate_ode_bdf;              <span class="co">// but use the stiff solver anyways</span></a>
<a class="sourceLine" id="cb9-43" title="43">  <span class="kw">using</span> stan::math::to_vector;</a>
<a class="sourceLine" id="cb9-44" title="44">  <span class="kw">using</span> stan::math::to_array_1d;</a>
<a class="sourceLine" id="cb9-45" title="45">  <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt; theta;</a>
<a class="sourceLine" id="cb9-46" title="46">  <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt; times = {<span class="dv">1</span>, <span class="dv">2</span>};</a>
<a class="sourceLine" id="cb9-47" title="47">  <span class="kw">auto</span> y = integrate_ode_bdf([&amp;A](<span class="kw">auto</span> t, <span class="kw">auto</span> y, </a>
<a class="sourceLine" id="cb9-48" title="48">                                  <span class="kw">auto</span> theta, <span class="kw">auto</span> x_r, <span class="kw">auto</span> x_i, <span class="bu">std::</span>ostream *msgs) {</a>
<a class="sourceLine" id="cb9-49" title="49">    <span class="cf">return</span> to_array_1d( (A * to_vector(y)).eval() );</a>
<a class="sourceLine" id="cb9-50" title="50">  }, to_array_1d(y0), <span class="dv">0</span>, times, theta, {}, {});</a>
<a class="sourceLine" id="cb9-51" title="51">  Eigen::VectorXd truth = stan::math::matrix_exp(A) * y0;</a>
<a class="sourceLine" id="cb9-52" title="52">  <span class="cf">return</span> (to_vector(y[<span class="dv">0</span>]) - truth).eval(); <span class="co">// should be &quot;zero&quot;</span></a>
<a class="sourceLine" id="cb9-53" title="53">}</a></code></pre></div>
<p>Again, in this compiled RMarkdown document, the <strong>knitr</strong> package has exported the <code>Cauchy</code>, <code>nonstiff</code> and <code>stiff</code> functions to R’s global environment using the <code>sourceCpp</code> function in the <strong>Rcpp</strong> package so that they can be called from R.</p>
<p>First, we numerically integrate the Cauchy PDF over its interquartile range — which has an area of <span class="math inline">\(\frac{1}{2}\)</span> — that we then double to verify that it is almost within machine precision of <span class="math inline">\(1\)</span>.</p>
<div class="sourceCode" id="cb10"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb10-1" title="1"><span class="kw">all.equal</span>(<span class="dv">1</span>, <span class="kw">Cauchy</span>(<span class="kw">rexp</span>(<span class="dv">1</span>)), <span class="dt">tol =</span> <span class="fl">1e-15</span>)</a>
<a class="sourceLine" id="cb10-2" title="2"><span class="co">#&gt; [1] TRUE</span></a></code></pre></div>
<p>Next, we consider the system of differential equations <span class="math display">\[\frac{d}{dt}\mathbf{y} = \mathbf{A}\mathbf{y}\]</span> where <span class="math inline">\(\mathbf{A}\)</span> is a square matrix such as that for a simple harmonic oscillator</p>
<p><span class="math display">\[\mathbf{A} = \begin{bmatrix}0 &amp; 1 \\ -1 &amp; -\theta\end{bmatrix}\]</span> for <span class="math inline">\(\theta \in \left(0,1\right)\)</span>. The solution for <span class="math inline">\(\mathbf{y}_t = e^{t\mathbf{A}}\mathbf{y}_0\)</span> can be obtained via the matrix exponential function, which is available in the Stan Math Library, but it can also be obtained numerically using a fourth-order Runge-Kutta solver, which is appropriate for non-stiff systems of ODEs, such as this one. However, it is possible, albeit less efficient in this case, to use the backward-differentiation formula solver for stiff systems of ODEs. In both cases, we calculate the difference between the analytical solution and the numerical one, and the stiff version does produce somewhat better accuracy in this case.</p>
<div class="sourceCode" id="cb11"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb11-1" title="1">A &lt;-<span class="st"> </span><span class="kw">matrix</span>(<span class="kw">c</span>(<span class="dv">0</span>, <span class="dv">-1</span>, <span class="dv">1</span>, <span class="op">-</span><span class="kw">runif</span>(<span class="dv">1</span>)), <span class="dt">nrow =</span> <span class="dv">2</span>, <span class="dt">ncol =</span> <span class="dv">2</span>)</a>
<a class="sourceLine" id="cb11-2" title="2">y0 &lt;-<span class="st"> </span><span class="kw">rexp</span>(<span class="dv">2</span>)</a>
<a class="sourceLine" id="cb11-3" title="3"><span class="kw">all.equal</span>(<span class="kw">nonstiff</span>(A, y0), <span class="kw">c</span>(<span class="dv">0</span>, <span class="dv">0</span>), <span class="dt">tol =</span> <span class="fl">1e-5</span>)</a>
<a class="sourceLine" id="cb11-4" title="4"><span class="co">#&gt; [1] TRUE</span></a>
<a class="sourceLine" id="cb11-5" title="5"><span class="kw">all.equal</span>(   <span class="kw">stiff</span>(A, y0), <span class="kw">c</span>(<span class="dv">0</span>, <span class="dv">0</span>), <span class="dt">tol =</span> <span class="fl">1e-8</span>)</a>
<a class="sourceLine" id="cb11-6" title="6"><span class="co">#&gt; [1] TRUE</span></a></code></pre></div>
</div>
<div id="parellelization" class="section level1">
<h1>Parellelization</h1>
<div id="map_rect-function" class="section level2">
<h2><code>map_rect</code> Function</h2>
<p>The Stan Math Library includes the <code>map_rect</code> function, which applies a function to each element of rectangular arrays and returns a vector, making it a bit like a restricted version of R’s <code>sapply</code> function. However, <code>map_rect</code> can also be executed in parallel by defining the pre-processor directive <code>STAN_THREADS</code> and then setting the <code>STAN_NUM_THREADS</code> environmental variable to be the number of threads to use, as in</p>
<div class="sourceCode" id="cb12"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb12-1" title="1"><span class="kw">Sys.setenv</span>(<span class="dt">STAN_NUM_THREADS =</span> <span class="dv">2</span>) <span class="co"># specify -1 to use all available cores</span></a></code></pre></div>
<p>Below is C++ code to test whether an integer is prime, using a rather brute-force algorithm and running it in parallel via <code>map_rect</code>.</p>
<div class="sourceCode" id="cb13"><pre class="sourceCode cpp"><code class="sourceCode cpp"><a class="sourceLine" id="cb13-1" title="1"><span class="co">// [[Rcpp::depends(BH)]]</span></a>
<a class="sourceLine" id="cb13-2" title="2"><span class="co">// [[Rcpp::depends(RcppEigen)]]</span></a>
<a class="sourceLine" id="cb13-3" title="3"><span class="co">// [[Rcpp::depends(RcppParallel)]]</span></a>
<a class="sourceLine" id="cb13-4" title="4"><span class="co">// [[Rcpp::depends(StanHeaders)]]</span></a>
<a class="sourceLine" id="cb13-5" title="5"><span class="pp">#include </span><span class="im">&lt;stan/math.hpp&gt;</span><span class="pp">                         </span><span class="co">// pulls in everything from rev/ and prim/</span></a>
<a class="sourceLine" id="cb13-6" title="6"><span class="pp">#include </span><span class="im">&lt;Rcpp.h&gt;</span></a>
<a class="sourceLine" id="cb13-7" title="7"><span class="pp">#include </span><span class="im">&lt;RcppEigen.h&gt;</span><span class="pp">                           </span><span class="co">// do this AFTER including stan/math</span></a>
<a class="sourceLine" id="cb13-8" title="8"></a>
<a class="sourceLine" id="cb13-9" title="9"><span class="co">// [[Rcpp::plugins(cpp17)]]</span></a>
<a class="sourceLine" id="cb13-10" title="10"></a>
<a class="sourceLine" id="cb13-11" title="11"><span class="co">// see https://en.wikipedia.org/wiki/Primality_test#Pseudocode</span></a>
<a class="sourceLine" id="cb13-12" title="12"><span class="kw">struct</span> is_prime {</a>
<a class="sourceLine" id="cb13-13" title="13">  is_prime() {}</a>
<a class="sourceLine" id="cb13-14" title="14">  <span class="kw">template</span> &lt;<span class="kw">typename</span> T1, <span class="kw">typename</span> T2&gt;</a>
<a class="sourceLine" id="cb13-15" title="15">  <span class="kw">auto</span></a>
<a class="sourceLine" id="cb13-16" title="16">  <span class="kw">operator</span>()(<span class="at">const</span> Eigen::Matrix&lt;T1, Eigen::Dynamic, <span class="dv">1</span>&gt;&amp; eta,</a>
<a class="sourceLine" id="cb13-17" title="17">             <span class="at">const</span> Eigen::Matrix&lt;T2, Eigen::Dynamic, <span class="dv">1</span>&gt;&amp; theta,</a>
<a class="sourceLine" id="cb13-18" title="18">             <span class="at">const</span> <span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt;&amp; x_r, <span class="at">const</span> <span class="bu">std::</span>vector&lt;<span class="dt">int</span>&gt;&amp; x_i,</a>
<a class="sourceLine" id="cb13-19" title="19">             <span class="bu">std::</span>ostream* msgs = <span class="kw">nullptr</span>) <span class="at">const</span> {</a>
<a class="sourceLine" id="cb13-20" title="20">    Eigen::VectorXd res(<span class="dv">1</span>); <span class="co">// can only return double or var vectors</span></a>
<a class="sourceLine" id="cb13-21" title="21">    <span class="dt">int</span> n = x_i[<span class="dv">0</span>];</a>
<a class="sourceLine" id="cb13-22" title="22">    <span class="cf">if</span> (n &lt;= <span class="dv">3</span>) {</a>
<a class="sourceLine" id="cb13-23" title="23">      res.coeffRef(<span class="dv">0</span>) = n &gt; <span class="dv">1</span>;</a>
<a class="sourceLine" id="cb13-24" title="24">      <span class="cf">return</span> res;</a>
<a class="sourceLine" id="cb13-25" title="25">    } <span class="cf">else</span> <span class="cf">if</span> ( (n % <span class="dv">2</span> == <span class="dv">0</span>) || (n % <span class="dv">3</span> == <span class="dv">0</span>) ) {</a>
<a class="sourceLine" id="cb13-26" title="26">      res.coeffRef(<span class="dv">0</span>) = <span class="kw">false</span>;</a>
<a class="sourceLine" id="cb13-27" title="27">      <span class="cf">return</span> res;</a>
<a class="sourceLine" id="cb13-28" title="28">    }</a>
<a class="sourceLine" id="cb13-29" title="29">    <span class="dt">int</span> i = <span class="dv">5</span>;</a>
<a class="sourceLine" id="cb13-30" title="30">    <span class="cf">while</span> (i * i &lt;= n) {</a>
<a class="sourceLine" id="cb13-31" title="31">      <span class="cf">if</span> ( (n % i == <span class="dv">0</span>) || (n % (i + <span class="dv">2</span>) == <span class="dv">0</span>) ) {</a>
<a class="sourceLine" id="cb13-32" title="32">        res.coeffRef(<span class="dv">0</span>) = <span class="kw">false</span>;</a>
<a class="sourceLine" id="cb13-33" title="33">        <span class="cf">return</span> res;</a>
<a class="sourceLine" id="cb13-34" title="34">      }</a>
<a class="sourceLine" id="cb13-35" title="35">      i += <span class="dv">6</span>;</a>
<a class="sourceLine" id="cb13-36" title="36">    }</a>
<a class="sourceLine" id="cb13-37" title="37">    res.coeffRef(<span class="dv">0</span>) = <span class="kw">true</span>;</a>
<a class="sourceLine" id="cb13-38" title="38">    <span class="cf">return</span> res;</a>
<a class="sourceLine" id="cb13-39" title="39">  }</a>
<a class="sourceLine" id="cb13-40" title="40">};</a>
<a class="sourceLine" id="cb13-41" title="41"></a>
<a class="sourceLine" id="cb13-42" title="42"><span class="co">/* parallelization */</span></a>
<a class="sourceLine" id="cb13-43" title="43"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb13-44" title="44"><span class="kw">auto</span> psapply(<span class="bu">std::</span>vector&lt;<span class="bu">std::</span>vector&lt;<span class="dt">int</span>&gt; &gt; n) {</a>
<a class="sourceLine" id="cb13-45" title="45">  <span class="bu">std::</span>vector&lt;Eigen::VectorXd&gt; eta(n.size()); <span class="co">// these all have to be the same size</span></a>
<a class="sourceLine" id="cb13-46" title="46">  Eigen::VectorXd theta;</a>
<a class="sourceLine" id="cb13-47" title="47">  <span class="bu">std::</span>vector&lt;<span class="bu">std::</span>vector&lt;<span class="dt">double</span>&gt; &gt; x_d(n.size());</a>
<a class="sourceLine" id="cb13-48" title="48">  <span class="cf">return</span> stan::math::map_rect&lt;<span class="dv">0</span>, is_prime&gt;(theta, eta, x_d, n, &amp;Rcpp::Rcout);</a>
<a class="sourceLine" id="cb13-49" title="49">}</a></code></pre></div>
<p>Since the signature for <code>n</code> is a <code>std::vector&lt;std::vector&lt;int&gt; &gt;</code>, we have to pass it from R as a list (which is converted to the outer <code>std::vector&lt;&gt;</code>) of integer vectors (which is converted to the inner <code>std::vector&lt;int&gt;</code>) that happen to be of size one in this case.</p>
<div class="sourceCode" id="cb14"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb14-1" title="1">odd &lt;-<span class="st"> </span><span class="kw">seq.int</span>(<span class="dt">from =</span> <span class="dv">2</span><span class="op">^</span><span class="dv">25</span> <span class="op">-</span><span class="st"> </span><span class="dv">1</span>, <span class="dt">to =</span> <span class="dv">2</span><span class="op">^</span><span class="dv">26</span> <span class="op">-</span><span class="st"> </span><span class="dv">1</span>, <span class="dt">by =</span> <span class="dv">2</span>)</a>
<a class="sourceLine" id="cb14-2" title="2"><span class="kw">tail</span>(<span class="kw">psapply</span>(<span class="dt">n =</span> <span class="kw">as.list</span>(odd))) <span class="op">==</span><span class="st"> </span><span class="dv">1</span> <span class="co"># check your process manager while this is running</span></a>
<a class="sourceLine" id="cb14-3" title="3"><span class="co">#&gt; [1] FALSE FALSE FALSE  TRUE FALSE FALSE</span></a></code></pre></div>
<p>Thus, <span class="math inline">\(2^{26} - 5 = 67,108,859\)</span> is a prime number.</p>
</div>
<div id="reduce_sum-function" class="section level2">
<h2><code>reduce_sum</code> Function</h2>
<p>The <code>reduce_sum</code> function can be used to sum a function of recursively-partitioned data in parallel. The Probability Mass Function (PMF) of the logarithmic distribution is</p>
<p><span class="math display">\[\Pr\left(x = k \mid p\right) = \frac{p^k}{-k\ln\left(1 - p\right)} \]</span> for a positive integer <span class="math inline">\(k\)</span> and <span class="math inline">\(0 &lt; p &lt; 1\)</span>. To verify that this PMF sums to <span class="math inline">\(1\)</span> over the natural numbers, we could accumulate it for a large finite number of terms, but would like to do so in parallel. To do so, we need to define a <code>partial_sum</code> function that adds the PMF for some subset of natural numbers and pass it to the <code>reduce_sum</code> function like</p>
<div class="sourceCode" id="cb15"><pre class="sourceCode cpp"><code class="sourceCode cpp"><a class="sourceLine" id="cb15-1" title="1"><span class="co">// [[Rcpp::depends(BH)]]</span></a>
<a class="sourceLine" id="cb15-2" title="2"><span class="co">// [[Rcpp::depends(RcppEigen)]]</span></a>
<a class="sourceLine" id="cb15-3" title="3"><span class="co">// [[Rcpp::depends(RcppParallel)]]</span></a>
<a class="sourceLine" id="cb15-4" title="4"><span class="co">// [[Rcpp::depends(StanHeaders)]]</span></a>
<a class="sourceLine" id="cb15-5" title="5"><span class="pp">#include </span><span class="im">&lt;stan/math.hpp&gt;</span><span class="pp">                         </span><span class="co">// pulls in everything from rev/ and prim/</span></a>
<a class="sourceLine" id="cb15-6" title="6"><span class="pp">#include </span><span class="im">&lt;Rcpp.h&gt;</span></a>
<a class="sourceLine" id="cb15-7" title="7"><span class="pp">#include </span><span class="im">&lt;RcppEigen.h&gt;</span><span class="pp">                           </span><span class="co">// do this AFTER including stan/math</span></a>
<a class="sourceLine" id="cb15-8" title="8"></a>
<a class="sourceLine" id="cb15-9" title="9"><span class="co">// [[Rcpp::plugins(cpp17)]]</span></a>
<a class="sourceLine" id="cb15-10" title="10"></a>
<a class="sourceLine" id="cb15-11" title="11"><span class="kw">template</span> &lt;<span class="kw">typename</span> T&gt;</a>
<a class="sourceLine" id="cb15-12" title="12"><span class="kw">struct</span> partial_sum {</a>
<a class="sourceLine" id="cb15-13" title="13">  partial_sum() {}</a>
<a class="sourceLine" id="cb15-14" title="14">  T <span class="kw">operator</span>()(<span class="at">const</span> <span class="bu">std::</span>vector&lt;<span class="dt">int</span>&gt;&amp; k_slice, </a>
<a class="sourceLine" id="cb15-15" title="15">               <span class="dt">int</span> start, <span class="dt">int</span> end, <span class="co">// these are ignored in this example</span></a>
<a class="sourceLine" id="cb15-16" title="16">               <span class="bu">std::</span>ostream* msgs, <span class="dt">double</span> p) {</a>
<a class="sourceLine" id="cb15-17" title="17">    <span class="dt">double</span> S = <span class="dv">0</span>;</a>
<a class="sourceLine" id="cb15-18" title="18">    <span class="cf">for</span> (<span class="dt">int</span> n = <span class="dv">0</span>; n &lt; k_slice.size(); n++) {</a>
<a class="sourceLine" id="cb15-19" title="19">      <span class="dt">int</span> k = k_slice[n];</a>
<a class="sourceLine" id="cb15-20" title="20">      S += stan::math::pow(p, k) / k;</a>
<a class="sourceLine" id="cb15-21" title="21">    }</a>
<a class="sourceLine" id="cb15-22" title="22">    <span class="cf">return</span> S;</a>
<a class="sourceLine" id="cb15-23" title="23">  }</a>
<a class="sourceLine" id="cb15-24" title="24">};</a>
<a class="sourceLine" id="cb15-25" title="25"></a>
<a class="sourceLine" id="cb15-26" title="26"><span class="co">// [[Rcpp::export]]</span></a>
<a class="sourceLine" id="cb15-27" title="27"><span class="dt">double</span> check_logarithmic_PMF(<span class="at">const</span> <span class="dt">double</span> p, <span class="at">const</span> <span class="dt">int</span> N = <span class="dv">1000</span>) {</a>
<a class="sourceLine" id="cb15-28" title="28">  <span class="kw">using</span> stan::math::reduce_sum;</a>
<a class="sourceLine" id="cb15-29" title="29">  <span class="bu">std::</span>vector&lt;<span class="dt">int</span>&gt; k(N);</a>
<a class="sourceLine" id="cb15-30" title="30">  <span class="bu">std::</span>iota(k.begin(), k.end(), <span class="dv">1</span>);</a>
<a class="sourceLine" id="cb15-31" title="31">  <span class="cf">return</span> reduce_sum&lt;partial_sum&lt;<span class="dt">double</span>&gt; &gt;(k, <span class="dv">1</span>, <span class="kw">nullptr</span>, p) / -<span class="bu">std::</span>log1p(-p);</a>
<a class="sourceLine" id="cb15-32" title="32">}</a></code></pre></div>
<p>The second argument passed to <code>reduce_sum</code> is the <code>grainsize</code>, which governs the size (and number) of the recursive subsets of <code>k</code> that are passed to the <code>partial_sum</code> function. A <code>grainsize</code> of <span class="math inline">\(1\)</span> indicates that the software will try to automatically tune the subset size for good performance, but that may be worse than choosing a <code>grainsize</code> by hand in a problem-specific fashion.</p>
<p>In any event, we can call the wrapper function <code>check_logarithmic_PMF</code> in R to verify that it returns <span class="math inline">\(1\)</span> to numerical tolerance:</p>
<div class="sourceCode" id="cb16"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb16-1" title="1"><span class="kw">stopifnot</span>(<span class="kw">all.equal</span>(<span class="dv">1</span>, <span class="kw">check_logarithmic_PMF</span>(<span class="dt">p =</span> <span class="dv">1</span> <span class="op">/</span><span class="st"> </span><span class="kw">sqrt</span>(<span class="dv">2</span>))))</a></code></pre></div>
</div>
</div>
<div id="defining-a-stan-model-in-c" class="section level1">
<h1>Defining a Stan Model in C++</h1>
<p>The Stan <em>language</em> does not have much support for sparse matrices for a variety of reasons. Essentially the only applicable function is <code>csr_matrix_times_vector</code>, which pre-multiplies a vector by a sparse matrix in compressed row storage by taking as arguments its number of rows, columns, non-zero values, column indices of non-zero values, and locations where the non-zero values start in each row. While the <code>csr_matrix_times_vector</code> function could be used to implement the example below, we illustrate how to use the sparse data structures in the <strong>Matrix</strong> and <strong>RcppEigen</strong> packages in a Stan model written in C++, which could easily be extended to more complicated models with sparse data structures.</p>
<p>Our C++ file for the log-likelihood of a linear model with a sparse design matrix reads as</p>
<pre><code>#include &lt;stan/model/model_header.hpp&gt;
#include &lt;Rcpp.h&gt;
#include &lt;RcppEigen.h&gt;

class sparselm_stan {

public: // these would ordinarily be private in the C++ code generated by Stan
  Eigen::Map&lt;Eigen::SparseMatrix&lt;double&gt; &gt; X;
  Eigen::VectorXd y;

  sparselm_stan(Eigen::Map&lt;Eigen::SparseMatrix&lt;double&gt; &gt; X, Eigen::VectorXd y) :
    X(X), y(y) {}

  template &lt;bool propto__ = false, bool jacobian__ = false, typename T__ = double&gt;
  // propto__ is usually true but causes log_prob() to return 0 when called from R
  // jacobian__ is usually true for MCMC but typically is false for optimization
  T__ log_prob(std::vector&lt;T__&gt;&amp; params_r__) const {
    using namespace stan::math;
    T__ lp__(0.0);
    accumulator&lt;T__&gt; lp_accum__;

    // set up model parameters
    std::vector&lt;int&gt; params_i__;
    stan::io::deserializer&lt;T__&gt; in__(params_r__, params_i__);
    auto beta = in__.template read&lt;Eigen::Matrix&lt;T__,-1,1&gt; &gt;(X.cols());
    auto sigma = in__.template read_constrain_lb&lt;T__, jacobian__&gt;(0, lp__);

    // log-likelihood (should add priors)
    lp_accum__.add(lp__);
    lp_accum__.add(normal_lpdf&lt;propto__&gt;(y, (X * beta).eval(), sigma));
    return lp_accum__.sum();
  }

  template &lt;bool propto__ = false, bool jacobian__ = false&gt;
  std::vector&lt;double&gt; gradient(std::vector&lt;double&gt;&amp; params_r__) const {
    // Calculate gradients using reverse-mode autodiff
    // although you could do them analytically in this case

    using std::vector;
    using stan::math::var;
    double lp;
    std::vector&lt;double&gt; gradient;
    try {
      vector&lt;var&gt; ad_params_r(params_r__.size());
      for (size_t i = 0; i &lt; params_r__.size(); ++i) {
        var var_i(params_r__[i]);
        ad_params_r[i] = var_i;
      }
      var adLogProb
        = this-&gt;log_prob&lt;propto__, jacobian__&gt;(ad_params_r);
      lp = adLogProb.val();
      adLogProb.grad(ad_params_r, gradient);
    } catch (const std::exception &amp;ex) {
      stan::math::recover_memory();
      throw;
    }
    stan::math::recover_memory();
    return gradient;
  }
};</code></pre>
<p>To use it from R, we call the <code>exposeClass</code> function in the <strong>Rcpp</strong> package with the necessary arguments and then call <code>sourceCpp</code> on the file it wrote in the temporary directory:</p>
<div class="sourceCode" id="cb18"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb18-1" title="1"><span class="kw">library</span>(Rcpp)</a>
<a class="sourceLine" id="cb18-2" title="2">tf &lt;-<span class="st"> </span><span class="kw">tempfile</span>(<span class="dt">fileext =</span> <span class="st">&quot;Module.cpp&quot;</span>)</a>
<a class="sourceLine" id="cb18-3" title="3"><span class="kw">exposeClass</span>(<span class="st">&quot;sparselm_stan&quot;</span>,</a>
<a class="sourceLine" id="cb18-4" title="4">      <span class="dt">constructors =</span> <span class="kw">list</span>(<span class="kw">c</span>(<span class="st">&quot;Eigen::Map&lt;Eigen::SparseMatrix&lt;double&gt; &gt;&quot;</span>, </a>
<a class="sourceLine" id="cb18-5" title="5">                            <span class="st">&quot;Eigen::VectorXd&quot;</span>)),</a>
<a class="sourceLine" id="cb18-6" title="6">      <span class="dt">fields =</span> <span class="kw">c</span>(<span class="st">&quot;X&quot;</span>, <span class="st">&quot;y&quot;</span>),</a>
<a class="sourceLine" id="cb18-7" title="7">      <span class="dt">methods =</span> <span class="kw">c</span>(<span class="st">&quot;log_prob&lt;&gt;&quot;</span>, <span class="st">&quot;gradient&lt;&gt;&quot;</span>),</a>
<a class="sourceLine" id="cb18-8" title="8">      <span class="dt">rename =</span> <span class="kw">c</span>(<span class="dt">log_prob =</span> <span class="st">&quot;log_prob&lt;&gt;&quot;</span>, <span class="dt">gradient =</span> <span class="st">&quot;gradient&lt;&gt;&quot;</span>),</a>
<a class="sourceLine" id="cb18-9" title="9">      <span class="dt">header =</span> <span class="kw">c</span>(<span class="st">&quot;// [[Rcpp::depends(BH)]]&quot;</span>,</a>
<a class="sourceLine" id="cb18-10" title="10">                 <span class="st">&quot;// [[Rcpp::depends(RcppEigen)]]&quot;</span>,</a>
<a class="sourceLine" id="cb18-11" title="11">                 <span class="st">&quot;// [[Rcpp::depends(RcppParallel)]&quot;</span>,</a>
<a class="sourceLine" id="cb18-12" title="12">                 <span class="st">&quot;// [[Rcpp::depends(StanHeaders)]]&quot;</span>,</a>
<a class="sourceLine" id="cb18-13" title="13">                 <span class="st">&quot;// [[Rcpp::plugins(cpp17)]]&quot;</span>,</a>
<a class="sourceLine" id="cb18-14" title="14">                 <span class="kw">paste0</span>(<span class="st">&quot;#include &lt;&quot;</span>, <span class="kw">file.path</span>(<span class="kw">getwd</span>(), <span class="st">&quot;sparselm_stan.hpp&quot;</span>), <span class="st">&quot;&gt;&quot;</span>)),</a>
<a class="sourceLine" id="cb18-15" title="15">      <span class="dt">file =</span> tf,</a>
<a class="sourceLine" id="cb18-16" title="16">      <span class="dt">Rfile =</span> <span class="ot">FALSE</span>)</a>
<a class="sourceLine" id="cb18-17" title="17"><span class="kw">Sys.setenv</span>(<span class="dt">PKG_CXXFLAGS =</span> <span class="kw">paste0</span>(<span class="kw">Sys.getenv</span>(<span class="st">&quot;PKG_CXXFLAGS&quot;</span>), <span class="st">&quot; -I&quot;</span>,</a>
<a class="sourceLine" id="cb18-18" title="18">                                 <span class="kw">system.file</span>(<span class="st">&quot;include&quot;</span>, <span class="st">&quot;src&quot;</span>, </a>
<a class="sourceLine" id="cb18-19" title="19">                                             <span class="dt">package =</span> <span class="st">&quot;StanHeaders&quot;</span>, <span class="dt">mustWork =</span> <span class="ot">TRUE</span>)))</a>
<a class="sourceLine" id="cb18-20" title="20"><span class="kw">sourceCpp</span>(tf)</a>
<a class="sourceLine" id="cb18-21" title="21">sparselm_stan</a>
<a class="sourceLine" id="cb18-22" title="22"><span class="co">#&gt; C++ class 'sparselm_stan' &lt;0x557ae52e0380&gt;</span></a>
<a class="sourceLine" id="cb18-23" title="23"><span class="co">#&gt; Constructors:</span></a>
<a class="sourceLine" id="cb18-24" title="24"><span class="co">#&gt;     sparselm_stan(Eigen::Map&lt;Eigen::SparseMatrix&lt;double, 0, int&gt;, 0, Eigen::Stride&lt;0, 0&gt; &gt;, Eigen::Matrix&lt;double, -1, 1, 0, -1, 1&gt;)</span></a>
<a class="sourceLine" id="cb18-25" title="25"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb18-26" title="26"><span class="co">#&gt; Fields: </span></a>
<a class="sourceLine" id="cb18-27" title="27"><span class="co">#&gt;     Eigen::Map&lt;Eigen::SparseMatrix&lt;double, 0, int&gt;, 0, Eigen::Stride&lt;0, 0&gt; &gt; X</span></a>
<a class="sourceLine" id="cb18-28" title="28"><span class="co">#&gt;     Eigen::Matrix&lt;double, -1, 1, 0, -1, 1&gt; y</span></a>
<a class="sourceLine" id="cb18-29" title="29"><span class="co">#&gt; </span></a>
<a class="sourceLine" id="cb18-30" title="30"><span class="co">#&gt; Methods: </span></a>
<a class="sourceLine" id="cb18-31" title="31"><span class="co">#&gt;      std::vector&lt;double, std::allocator&lt;double&gt; &gt; gradient(std::vector&lt;double, std::allocator&lt;double&gt; &gt;)  const </span></a>
<a class="sourceLine" id="cb18-32" title="32"><span class="co">#&gt;            </span></a>
<a class="sourceLine" id="cb18-33" title="33"><span class="co">#&gt;      double log_prob(std::vector&lt;double, std::allocator&lt;double&gt; &gt;)  const </span></a>
<a class="sourceLine" id="cb18-34" title="34"><span class="co">#&gt; </span></a></code></pre></div>
<p>At this point, we need a sparse design matrix and (dense) outcome vector to pass to the constructor. The former can be created with a variety of functions in the <strong>Matrix</strong> package, such as</p>
<div class="sourceCode" id="cb19"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb19-1" title="1">dd &lt;-<span class="st"> </span><span class="kw">data.frame</span>(<span class="dt">a =</span> <span class="kw">gl</span>(<span class="dv">3</span>, <span class="dv">4</span>), <span class="dt">b =</span> <span class="kw">gl</span>(<span class="dv">4</span>, <span class="dv">1</span>, <span class="dv">12</span>))</a>
<a class="sourceLine" id="cb19-2" title="2">X &lt;-<span class="st"> </span>Matrix<span class="op">::</span><span class="kw">sparse.model.matrix</span>(<span class="op">~</span><span class="st"> </span>a <span class="op">+</span><span class="st"> </span>b, <span class="dt">data =</span> dd)</a>
<a class="sourceLine" id="cb19-3" title="3">X</a>
<a class="sourceLine" id="cb19-4" title="4"><span class="co">#&gt; 12 x 6 sparse Matrix of class &quot;dgCMatrix&quot;</span></a>
<a class="sourceLine" id="cb19-5" title="5"><span class="co">#&gt;    (Intercept) a2 a3 b2 b3 b4</span></a>
<a class="sourceLine" id="cb19-6" title="6"><span class="co">#&gt; 1            1  .  .  .  .  .</span></a>
<a class="sourceLine" id="cb19-7" title="7"><span class="co">#&gt; 2            1  .  .  1  .  .</span></a>
<a class="sourceLine" id="cb19-8" title="8"><span class="co">#&gt; 3            1  .  .  .  1  .</span></a>
<a class="sourceLine" id="cb19-9" title="9"><span class="co">#&gt; 4            1  .  .  .  .  1</span></a>
<a class="sourceLine" id="cb19-10" title="10"><span class="co">#&gt; 5            1  1  .  .  .  .</span></a>
<a class="sourceLine" id="cb19-11" title="11"><span class="co">#&gt; 6            1  1  .  1  .  .</span></a>
<a class="sourceLine" id="cb19-12" title="12"><span class="co">#&gt; 7            1  1  .  .  1  .</span></a>
<a class="sourceLine" id="cb19-13" title="13"><span class="co">#&gt; 8            1  1  .  .  .  1</span></a>
<a class="sourceLine" id="cb19-14" title="14"><span class="co">#&gt; 9            1  .  1  .  .  .</span></a>
<a class="sourceLine" id="cb19-15" title="15"><span class="co">#&gt; 10           1  .  1  1  .  .</span></a>
<a class="sourceLine" id="cb19-16" title="16"><span class="co">#&gt; 11           1  .  1  .  1  .</span></a>
<a class="sourceLine" id="cb19-17" title="17"><span class="co">#&gt; 12           1  .  1  .  .  1</span></a></code></pre></div>
<p>Finally, we call the <code>new</code> function in the <strong>methods</strong> package, which essentially calls our C++ constructor and provides an R interface to the instantiated object, which contains the <code>log_prob</code> and <code>gradient</code> methods we defined and can be called with arbitrary inputs.</p>
<div class="sourceCode" id="cb20"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb20-1" title="1">sm &lt;-<span class="st"> </span><span class="kw">new</span>(sparselm_stan, <span class="dt">X =</span> X, <span class="dt">y =</span> <span class="kw">rnorm</span>(<span class="kw">nrow</span>(X)))</a>
<a class="sourceLine" id="cb20-2" title="2">sm<span class="op">$</span><span class="kw">log_prob</span>(<span class="kw">c</span>(<span class="dt">beta =</span> <span class="kw">rnorm</span>(<span class="kw">ncol</span>(X)), <span class="dt">log_sigma =</span> <span class="kw">log</span>(pi)))</a>
<a class="sourceLine" id="cb20-3" title="3"><span class="co">#&gt; [1] -26.01577</span></a>
<a class="sourceLine" id="cb20-4" title="4"><span class="kw">round</span>(sm<span class="op">$</span><span class="kw">gradient</span>(<span class="kw">c</span>(<span class="dt">beta =</span> <span class="kw">rnorm</span>(<span class="kw">ncol</span>(X)), <span class="dt">log_sigma =</span> <span class="kw">log</span>(pi))), <span class="dt">digits =</span> <span class="dv">4</span>)</a>
<a class="sourceLine" id="cb20-5" title="5"><span class="co">#&gt; [1] -2.1984 -1.0424 -0.5089 -0.4552 -0.6047 -0.8303 -6.6551</span></a></code></pre></div>
</div>



<!-- code folding -->


<!-- dynamically load mathjax for compatibility with self-contained -->
<script>
  (function () {
    var script = document.createElement("script");
    script.type = "text/javascript";
    script.src  = "https://mathjax.rstudio.com/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML";
    document.getElementsByTagName("head")[0].appendChild(script);
  })();
</script>

</body>
</html>