<|start_filename|>docs/reference/tsne_plot.html<|end_filename|>
<!-- Generated by pkgdown: do not edit by hand -->
<!DOCTYPE html>
<html lang="en">
  <head>
  <meta charset="utf-8">
<meta http-equiv="X-UA-Compatible" content="IE=edge">
<meta name="viewport" content="width=device-width, initial-scale=1.0">

<title>t-SNE Plot — tsne_plot • fastTopics</title>


<!-- jquery -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/jquery/3.4.1/jquery.min.js" integrity="sha256-CSXorXvZcTkaix6Yvo6HppcZGetbYMGWSFlBw8HfCJo=" crossorigin="anonymous"></script>
<!-- Bootstrap -->

<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/css/bootstrap.min.css" integrity="sha256-bZLfwXAP04zRMK2BjiO8iu9pf4FbLqX6zitd+tIvLhE=" crossorigin="anonymous" />

<script src="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/js/bootstrap.min.js" integrity="sha256-nuL8/2cJ5NDSSwnKD8VqreErSWHtnEP9E7AySL+1ev4=" crossorigin="anonymous"></script>

<!-- bootstrap-toc -->
<link rel="stylesheet" href="../bootstrap-toc.css">
<script src="../bootstrap-toc.js"></script>

<!-- Font Awesome icons -->
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/all.min.css" integrity="sha256-mmgLkCYLUQbXn0B1SRqzHar6dCnv9oZFPEC1g1cwlkk=" crossorigin="anonymous" />
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/v4-shims.min.css" integrity="sha256-wZjR52fzng1pJHwx4aV2AO3yyTOXrcDW7jBpJtTwVxw=" crossorigin="anonymous" />

<!-- clipboard.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/clipboard.js/2.0.6/clipboard.min.js" integrity="sha256-inc5kl9MA1hkeYUt+EC3BhlIgyp/2jDIyBLS6k3UxPI=" crossorigin="anonymous"></script>

<!-- headroom.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/headroom.min.js" integrity="sha256-AsUX4SJE1+yuDu5+mAVzJbuYNPHj/WroHuZ8Ir/CkE0=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/jQuery.headroom.min.js" integrity="sha256-ZX/yNShbjqsohH1k95liqY9Gd8uOiE1S4vZc+9KQ1K4=" crossorigin="anonymous"></script>

<!-- pkgdown -->
<link href="../pkgdown.css" rel="stylesheet">
<script src="../pkgdown.js"></script>




<meta property="og:title" content="t-SNE Plot — tsne_plot" />
<meta property="og:description" content="Visualize the &quot;structure&quot; of the Poisson NMF loadings
  or multinomial topic model mixture proportions by projection onto a
  2-d surface. Samples in the projection are colored according to the
  their loadings/proportions. By default, t-SNE is used to compute
  the 2-d embedding from the loadings or mixture proportions." />




<!-- mathjax -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/MathJax.js" integrity="sha256-nvJJv9wWKEm88qvoQl9ekL2J+k/RWIsaSScxxlsrv8k=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/config/TeX-AMS-MML_HTMLorMML.js" integrity="sha256-84DKXVJXs0/F8OTMzX4UR909+jtl4G7SPypPavF+GfA=" crossorigin="anonymous"></script>

<!--[if lt IE 9]>
<script src="https://oss.maxcdn.com/html5shiv/3.7.3/html5shiv.min.js"></script>
<script src="https://oss.maxcdn.com/respond/1.4.2/respond.min.js"></script>
<![endif]-->



  </head>

  <body data-spy="scroll" data-target="#toc">
    <div class="container template-reference-topic">
      <header>
      <div class="navbar navbar-default navbar-fixed-top" role="navigation">
  <div class="container">
    <div class="navbar-header">
      <button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#navbar" aria-expanded="false">
        <span class="sr-only">Toggle navigation</span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
      </button>
      <span class="navbar-brand">
        <a class="navbar-link" href="../index.html">fastTopics</a>
        <span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.5-15</span>
      </span>
    </div>

    <div id="navbar" class="navbar-collapse collapse">
      <ul class="nav navbar-nav">
        <li>
  <a href="../index.html">Home</a>
</li>
<li>
  <a href="../articles/index.html">Vignettes</a>
</li>
<li>
  <a href="../reference/index.html">Functions</a>
</li>
      </ul>
      <ul class="nav navbar-nav navbar-right">
        <li>
  <a href="https://github.com/stephenslab/fastTopics">Source</a>
</li>
      </ul>
      
    </div><!--/.nav-collapse -->
  </div><!--/.container -->
</div><!--/.navbar -->

      

      </header>

<div class="row">
  <div class="col-md-9 contents">
    <div class="page-header">
    <h1>t-SNE Plot</h1>
    <small class="dont-index">Source: <a href='https://github.com/stephenslab/fastTopics/blob/master/R/plots.R'><code>R/plots.R</code></a></small>
    <div class="hidden name"><code>tsne_plot.Rd</code></div>
    </div>

    <div class="ref-description">
    <p>Visualize the "structure" of the Poisson NMF loadings
  or multinomial topic model mixture proportions by projection onto a
  2-d surface. Samples in the projection are colored according to the
  their loadings/proportions. By default, t-SNE is used to compute
  the 2-d embedding from the loadings or mixture proportions.</p>
    </div>

    <pre class="usage"><span class='fu'>tsne_plot</span>(
  <span class='no'>fit</span>,
  <span class='kw'>color</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/c.html'>c</a></span>(<span class='st'>"mixprop"</span>, <span class='st'>"loading"</span>),
  <span class='no'>k</span>,
  <span class='no'>tsne</span>,
  <span class='kw'>ggplot_call</span> <span class='kw'>=</span> <span class='no'>tsne_plot_ggplot_call</span>,
  <span class='kw'>plot_grid_call</span> <span class='kw'>=</span> <span class='kw'>function</span>(<span class='no'>plots</span>) <span class='fu'><a href='https://rdrr.io/r/base/do.call.html'>do.call</a></span>(<span class='no'>plot_grid</span>, <span class='no'>plots</span>),
  <span class='no'>...</span>
)

<span class='fu'>tsne_plot_ggplot_call</span>(<span class='no'>dat</span>, <span class='no'>topic.label</span>, <span class='kw'>font.size</span> <span class='kw'>=</span> <span class='fl'>9</span>)</pre>

    <h2 class="hasAnchor" id="arguments"><a class="anchor" href="#arguments"></a>Arguments</h2>
    <table class="ref-arguments">
    <colgroup><col class="name" /><col class="desc" /></colgroup>
    <tr>
      <th>fit</th>
      <td><p>An object of class &#8220;poisson_nmf_fit&#8221; or
&#8220;multinom_topic_model_fit&#8221;.</p></td>
    </tr>
    <tr>
      <th>color</th>
      <td><p>The data mapped to the color &#8220;aesthetic&#8221; in the
plot. When <code>color = "loading"</code>, the estimated loadings (stored
the <code>L</code> matrix) in the Poisson NMF model are plotted; when
<code>color = "mixprop"</code>, the estimated mixture proportions (which are
recovered from the loadings by calling
<code><a href='poisson2multinom.html'>poisson2multinom</a></code> are shown. When <code>fit</code> is a
&#8220;multinom_topic_model_fit&#8221; object, the only available option
is <code>color = "mixprop"</code>. In most settings, the mixture proportions
are preferred, even if the 2-d embedding is computed from the
loading matrix of the Poisson NMF model.</p></td>
    </tr>
    <tr>
      <th>k</th>
      <td><p>The topic, or topics, selected by number or name. One plot
is created per selected topic. When not specified, all topics are
plotted.</p></td>
    </tr>
    <tr>
      <th>tsne</th>
      <td><p>A 2-d embedding of the samples (rows of X), or a subset
of the samples, such as an output from
<code><a href='tsne_from_topics.html'>tsne_from_topics</a></code>. It should be a list object with the
same structure as a <code>tsne_from_topics</code> output; see
<code><a href='tsne_from_topics.html'>tsne_from_topics</a></code> for details. If not provided, a 2-d
t-SNE embedding will be estimated automatically by calling
<code>tsne_from_topics</code>.</p></td>
    </tr>
    <tr>
      <th>ggplot_call</th>
      <td><p>The function used to create the plot. Replace
<code>tsne_plot_ggplot_call</code> with your own function to
customize the appearance of the plot.</p></td>
    </tr>
    <tr>
      <th>plot_grid_call</th>
      <td><p>When multiple topics are selected, this is
the function used to arrange the plots into a grid using
<code><a href='https://rdrr.io/pkg/cowplot/man/plot_grid.html'>plot_grid</a></code>. It should be a function accepting
a single argument, <code>plots</code>, a list of <code>ggplot</code> objects.</p></td>
    </tr>
    <tr>
      <th>...</th>
      <td><p>Additional arguments passed to
<code>tsne_from_topics</code>. These arguments are only used if
<code>tsne</code> is not provided.</p></td>
    </tr>
    <tr>
      <th>dat</th>
      <td><p>A data frame passed as input to
<code><a href='https://ggplot2.tidyverse.org/reference/ggplot.html'>ggplot</a></code>, containing, at a minimum, columns
&#8220;d1&#8221;, &#8220;d2&#8221; (the first and second dimensions in the
2-d embedding), and &#8220;loading&#8221;.</p></td>
    </tr>
    <tr>
      <th>topic.label</th>
      <td><p>The name or number of the topic being plotted;
it is only used to determine the plot title.</p></td>
    </tr>
    <tr>
      <th>font.size</th>
      <td><p>Font size used in plot.</p></td>
    </tr>
    </table>

    <h2 class="hasAnchor" id="value"><a class="anchor" href="#value"></a>Value</h2>

    <p>A <code>ggplot</code> object.</p>
    <h2 class="hasAnchor" id="details"><a class="anchor" href="#details"></a>Details</h2>

    <p>This is a lightweight interface primarily intended to
expedite creation of scatterplots for visualizing the loadings or
mixture proportions in 2-d; most of the &#8220;heavy lifting&#8221; is
done by ggplot2. The 2-d embedding itself is computed by invoking
function <code><a href='tsne_from_topics.html'>tsne_from_topics</a></code> (unless the &#8220;tsne&#8221;
input is provided). For more control over the plot's appearance,
the plot can be customized by modifying the <code>ggplot_call</code> and
<code>plot_grid_call</code> arguments.</p>
<p>An effective 2-d visualization may also necessitate some
fine-tunning of the t-SNE settings, such as the
&#8220;perplexity&#8221;, or the number of samples included in the
plot. The t-SNE settings can be controlled by the additional
arguments (...) passed to <code>tsne_from_topics</code>; see
<code><a href='tsne_from_topics.html'>tsne_from_topics</a></code> for details. Alternatively, a 2-d
embedding may be pre-computed, and passed as argument <code>tsne</code>
to <code>tsne_plot</code>.</p>

  </div>
  <div class="col-md-3 hidden-xs hidden-sm" id="pkgdown-sidebar">
    <nav id="toc" data-toggle="toc" class="sticky-top">
      <h2 data-toc-skip>Contents</h2>
    </nav>
  </div>
</div>


      <footer>
      <div class="copyright">
  <p>Developed by <NAME>, <NAME>, <NAME>, <NAME>.</p>
</div>

<div class="pkgdown">
  <p>Site built with <a href="https://pkgdown.r-lib.org/">pkgdown</a> 1.5.1.</p>
</div>

      </footer>
   </div>

  


  </body>
</html>




<|start_filename|>src/scd.cpp<|end_filename|>
#include <RcppParallel.h>
#include "misc.h"
#include "poismix.h"

using namespace arma;

// FUNCTION DECLARATIONS
// ---------------------
void scd_update_factors (const mat& A, const mat& W, mat& H, const vec& j,
			 unsigned int numiter, double e);

void scd_update_factors_sparse (const sp_mat& A, const mat& W, mat& H,
				const vec& j, unsigned int numiter, 
				double e);

// INLINE FUNCTION DEFINITIONS
// ---------------------------
// Perform one or more SCD updates for a single column of H, in which
// A is approximated by the matrix product W*H.
inline void scd_update_factor (const mat& A, const mat& W, mat& H,
			       unsigned int j, unsigned int numiter, 
			       double e) {
  H.col(j) = scd_kl_update(W,A.col(j),H.col(j),numiter,e);
}

// Perform one or more SCD updates for a single column of H, in which
// A is approximated by the matrix product W*H, and A is a sparse
// matrix. Here, "sumw" should contain the precomputed column sums of
// W; that is, sumw = colSums(W).
inline void scd_update_factor_sparse (const sp_mat& A, const mat& W,
				      const vec& sumw, mat& H, unsigned int j,
				      unsigned int numiter, double e) {
  vec          a = nonzeros(A.col(j));
  unsigned int n = a.n_elem;
  uvec         i(n);
  getcolnonzeros(A,i,j);
  H.col(j) = scd_kl_update(W.rows(i),sumw,a,H.col(j),numiter,e);
}

// CLASS DEFINITIONS
// -----------------
// This class is used to implement multithreaded computation of the
// loadings updates in scd_update_factors_parallel_rcpp.
struct scd_factor_updater : public RcppParallel::Worker {
  const mat&   A;
  const mat&   W;
  mat&         H;
  const vec&   j;
  unsigned int numiter;
  double       e;

  // This is used to create a scd_factor_updater object.
  scd_factor_updater (const mat& A, const mat& W, mat& H, const vec& j,
		      unsigned int numiter, double e) :
    A(A), W(W), H(H), j(j), numiter(numiter), e(e) { };

  // This function updates the factors for a given range of columns.
  void operator() (std::size_t begin, std::size_t end) {
    for (unsigned int i = begin; i < end; i++)
      scd_update_factor(A,W,H,j(i),numiter,e);
  }
};

// This class is used to implement multithreaded computation of the
// loadings updates in scd_update_factors_sparse_parallel_rcpp.
struct scd_factor_updater_sparse : public RcppParallel::Worker {
  const sp_mat& A;
  const mat&    W;
  vec           sumw;
  mat&          H;
  const vec&    j;
  unsigned int  numiter;
  double        e;

  // This is used to create a scd_factor_updater_sparse object.
  scd_factor_updater_sparse (const sp_mat& A, const mat& W, mat& H, 
			     const vec& j, unsigned int numiter, double e) :
    A(A), W(W), sumw(W.n_cols), H(H), j(j), numiter(numiter), e(e) {
    sumw = sum(W,0);
  };
 
  // This function updates the factors for a given range of columns.
  void operator() (std::size_t begin, std::size_t end) {
    for (unsigned int i = begin; i < end; i++)
      scd_update_factor_sparse(A,W,sumw,H,j(i),numiter,e);
  }
};

// FUNCTION DEFINITIONS
// --------------------
// Perform sequential co-ordinate descent (SCD) updates for the
// factors matrix (H).
//
// The inputs are: A, the m x n matrix to be decomposed, in which A is
// approximated as W*H; W, the initial estimate of the n x k loadings
// matrix; H, the k x m factors matrix; numiter, the number of
// inner-loop iterations to perform; and e, a non-negative scalar
// included in the computations to prevent NaNs due to division by
// zero.
//
// This implementation is adapted from the R and C++ code developed by
// <NAME> and <NAME>, which is available for download at
// https://github.com/linxihui/NNLM.
//
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat scd_update_factors_rcpp (const arma::mat& A, const arma::mat& W,
				   const arma::mat& H, const arma::vec& j,
				   unsigned int numiter, double e) {
  mat Hnew = H;
  scd_update_factors(A,W,Hnew,j,numiter,e);
  return Hnew;
}

// Perform sequential co-ordinate descent (SCD) updates for the
// factors matrix (H) when the data matrix (A) is sparse.
//
// [[Rcpp::export]]
arma::mat scd_update_factors_sparse_rcpp (const arma::sp_mat& A,
					  const arma::mat& W,
					  const arma::mat& H,
					  const arma::vec& j,
					  unsigned int numiter, 
					  double e) {
  mat Hnew = H;
  scd_update_factors_sparse(A,W,Hnew,j,numiter,e);
  return Hnew;
}

// This does the same thing as scd_update_factors_rcpp, except that
// Intel Threading Building Blocks (TBB) are used to update the
// loadings in parallel.
//
// [[Rcpp::depends(RcppParallel)]]
// [[Rcpp::export]]
arma::mat scd_update_factors_parallel_rcpp (const arma::mat& A, 
					    const arma::mat& W,
					    const arma::mat& H, 
					    const arma::vec& j,
					    unsigned int numiter, 
					    double e) {
  mat Hnew = H;
  scd_factor_updater worker(A,W,Hnew,j,numiter,e);
  parallelFor(0,j.n_elem,worker);
  return Hnew;
}

// This does the same thing as scd_update_factors_sparse_rcpp, except
// that Intel Threading Building Blocks (TBB) are used to update the
// loadings in parallel.
//
// [[Rcpp::export]]
arma::mat scd_update_factors_sparse_parallel_rcpp (const arma::sp_mat& A, 
						   const arma::mat& W,
						   const arma::mat& H, 
						   const arma::vec& j,
						   unsigned int numiter, 
						   double e) {
  mat Hnew = H;
  scd_factor_updater_sparse worker(A,W,Hnew,j,numiter,e);
  parallelFor(0,j.n_elem,worker);
  return Hnew;
}

// Iterate the SCD updates over all columns of H, in which A is
// approximated by the matrix product W*H.
void scd_update_factors (const mat& A, const mat& W, mat& H, const vec& j,
			 unsigned int numiter, double e) {
  unsigned int n = j.n_elem;
  for (unsigned int i = 0; i < n; i++)
    scd_update_factor(A,W,H,j(i),numiter,e);
}

// This is the same as scd_update_factors, except that the count data are
// stored as a sparse matrix.
void scd_update_factors_sparse (const sp_mat& A, const mat& W, mat& H,
				const vec& j, unsigned int numiter, double e) {
  unsigned int n    = j.n_elem;
  vec          sumw = sum(W,0);
  for (unsigned int i = 0; i < n; i++)
    scd_update_factor_sparse(A,W,sumw,H,j(i),numiter,e);
}


<|start_filename|>src/mixem.cpp<|end_filename|>
#include "mixem.h"
#include "misc.h"

using namespace arma;

// FUNCTION DECLARATIONS
// ---------------------
void mixem_update (const arma::mat& L1, const arma::vec& w,
		   arma::vec& x, arma::mat& P);

// FUNCTION DEFINITIONS
// --------------------
// This is mainly used for testing the mixem C++ function.
//
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::vec mixem_rcpp (const arma::mat& L, const arma::vec& w,
		      const arma::vec& x0, unsigned int numiter) {
  return mixem(L,w,x0,numiter);
}

// Compute a maximum-likelihood estimate (MLE) of the mixture
// proportions in the multinomial mixture model by iterating the EM
// updates for a fixed number of iterations.
//
// Input argument L is an n x m matrix with non-negative entries;
// input w is a vector of length n containing a non-negative "weight"
// associated with each row of L; input argument x0 is the initial
// estimate of the mixture proportions; input P is a matrix of the
// same dimension as L, and is used to store the posterior mixture
// assignment probabilities; and input "numiter" specifies the number
// of EM updates to perform.
//
// The return value is a vector of length m containing the updated
// mixture proportions.
//
// Note that x0 and L need not be normalized; they will automatically
// be normalized inside this function.
//
// Also note that it does not make sense to compute a MLE of the
// mixture proportions when n < 2 and/or when m < 2; mixem will supply
// a result in such cases, but the result will not be valid.
vec mixem (const mat& L, const vec& w, const vec& x0, unsigned int numiter) {
  mat L1 = L;
  mat P  = L;
  vec x  = x0;
  normalizecols(L1);
  mixem(L1,w,x,P,numiter);
  return x;
}

// Use this variant of mixem if you plan on using the same L matrix
// multiple times, or for calling mixem multiple times with matrices
// of the same dimension. In the first case, you can reuse the L1 and
// P matrices; in the latter case, you can reuse the P matrix.
//
// For the result to be valid, the matrix L1 should be normalized
// beforehand so that each column sums to 1. P should be a matrix of
// the same size as L1.
//
// Note that x need not be normalized; it will automatically be
// normalized inside this function.
//
// Also note that in this mixem variant, L1 and w do not need to
// contain all the data; any rows of L1 associated with zero weights
// have no effect, so only the vector of nonzero weights w, and the
// rows of L1 associated with those weights, need to be supplied.
void mixem (const mat& L1, const vec& w, vec& x, mat& P, 
            unsigned int numiter) {
  for (unsigned int i = 0; i < numiter; i++)
    mixem_update(L1,w,x,P);
}

// Perform a single EM update. For this update to be valid, the matrix
// L1 should be normalized beforehand so that each column sums to 1.
// Note that x need not be normalized; it will automatically be
// normalized inside this function.
void mixem_update (const mat& L1, const vec& w, vec& x, mat& P) {
  double e = 1e-15;

  // Normalize the "weights".
  vec w1 = w/sum(w);

  // Normalize the mixture proportions.
  x /= sum(x);
  
  // Compute the posterior mixture assignment probabilities. A small
  // number is added to the posterior probabilities to prevent any
  // divisions by zero. This is the "E step".
  P = L1;
  scalecols(P,x);
  normalizerowsbymax(P);
  P += e;
  normalizerows(P);
    
  // Update the mixture weights. This is the "M step".
  x = trans(P) * w1;
}

// Find the maximum-likelihood estimate (MLE) for the special case
// when only one of the weights (w) is positive. Here, L1 should be
// the column-normalized matrix, and i should be the index of the
// nonzero weight.
void mixture_one_nonzero (const mat& L1, unsigned int i, vec& x) {
  unsigned int j = index_max(L1.row(i));
  x.fill(0);
  x(j) = 1;
}


<|start_filename|>src/poismix.h<|end_filename|>
#ifndef INCLUDE_POISMIXEM
#define INCLUDE_POISMIXEM

#include <RcppArmadillo.h>

// FUNCTION DECLARATIONS
// ---------------------
arma::vec poismixem (const arma::mat& L, const arma::vec& w, 
		     const arma::vec& x0, unsigned int numiter);

void poismixem (const arma::mat& L1, const arma::vec& u, const arma::vec& w, 
		arma::vec& x, arma::mat& P, unsigned int numiter);

void poismixem (const arma::mat& L1, const arma::vec& u, const arma::vec& w,
		const arma::uvec& i, arma::vec& x, unsigned int numiter);

arma::vec scd_kl_update (const arma::mat& L, const arma::vec& w,
			 const arma::vec& x0, unsigned int numiter, 
			 double e);

arma::vec scd_kl_update (const arma::mat& L, const arma::vec& u,
			 const arma::vec& w, const arma::vec& x0,
			 unsigned int numiter, double e);

arma::vec ccd_kl_update (const arma::mat& L, const arma::vec& w,
			 const arma::vec& x0, double e);

arma::vec ccd_kl_update (const arma::mat& L, const arma::vec& u,
			 const arma::vec& w, const arma::vec& x0,
			 double e);

#endif


<|start_filename|>src/cost.h<|end_filename|>
#ifndef INCLUDE_COST
#define INCLUDE_COST

#include <RcppArmadillo.h>

// FUNCTION DECLARATIONS
// ---------------------
arma::vec cost        (const arma::mat& X, const arma::mat& A,
		       const arma::mat& B, double e, bool poisson);
arma::vec cost_sparse (const arma::sp_mat& X, const arma::mat& A,
	   	       const arma::mat& B, double e, bool poisson);

#endif


<|start_filename|>src/ccd.cpp<|end_filename|>
#include <RcppParallel.h>
#include "misc.h"
#include "poismix.h"

using namespace arma;

// FUNCTION DECLARATIONS
// ---------------------
void ccd_update_factors (const mat& V, const mat& W, mat& H, double e);

void ccd_update_factors_sparse (const sp_mat& V, const mat& W, mat& H, 
				double e);

// INLINE FUNCTION DEFINITIONS
// ---------------------------
// Perform a CCD update for a single column of H, in which V is
// approximated by the matrix product W*H.
inline void ccd_update_factor (const mat& V, const mat& W, mat& H,
			       unsigned int j, double e) {
  H.col(j) = ccd_kl_update(W,V.col(j),H.col(j),e);
}

// Perform a CCD update for a single column of H, in which V is
// approximated by the matrix product W*H, and V is a sparse matrix.
// Here, "sumw" should contain the precomputed column sums of W; that
// is, sumw = colSums(W).
inline void ccd_update_factor_sparse (const sp_mat& V, const mat& W,
				      const vec& sumw, mat& H, 
				      unsigned int j, double e) {
  vec          v = nonzeros(V.col(j));
  unsigned int n = v.n_elem;
  uvec         i(n);
  getcolnonzeros(V,i,j);
  H.col(j) = ccd_kl_update(W.rows(i),sumw,v,H.col(j),e);
}

// CLASS DEFINITIONS
// -----------------
// This class is used to implement multithreaded computation of the
// loadings updates in ccd_update_factors_parallel_rcpp.
struct ccd_factor_updater : public RcppParallel::Worker {
  const mat& V;
  const mat& W;
  mat&       H;
  double     e;

  // This is used to create a ccd_factor_updater object.
  ccd_factor_updater (const mat& V, const mat& W, mat& H, double e) :
    V(V), W(W), H(H), e(e) { };

  // This function updates the factors for a given range of columns.
  void operator() (std::size_t begin, std::size_t end) {
    for (unsigned int j = begin; j < end; j++)
      ccd_update_factor(V,W,H,j,e);
  }
};

// This class is used to implement multithreaded computation of the
// loadings updates in ccd_update_factors_sparse_parallel_rcpp.
struct ccd_factor_updater_sparse : public RcppParallel::Worker {
  const sp_mat& V;
  const mat&    W;
  vec           sumw;
  mat&          H;
  double        e;

  // This is used to create a ccd_factor_updater_sparse object.
  ccd_factor_updater_sparse (const sp_mat& V, const mat& W, mat& H, double e) :
    V(V), W(W), sumw(W.n_cols), H(H), e(e) {
    sumw = sum(W,0);
  };

  // This function updates the factors for a given range of columns.
  void operator() (std::size_t begin, std::size_t end) {
    for (unsigned int j = begin; j < end; j++)
      ccd_update_factor_sparse(V,W,sumw,H,j,e);
  }
};

// FUNCTION DEFINITIONS
// --------------------
// Perform a single cyclic co-ordinate descent (CCD) update for the
// factors matrix (H). 
//
// The inputs are: V, the n x m data matrix; W, the n x k loadings
// matrix; H, the initial estimate of the k x m factors matrix; and
// e, a non-negative scalar specifying the minimum value of the
// updated parameters.
//
// This implementation is adapted from the C++ code developed by
// <NAME> and <NAME>, which is available for download
// at www.cs.utexas.edu/~cjhsieh/nmf.
//
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat ccd_update_factors_rcpp (const arma::mat& V, const arma::mat& W, 
				   arma::mat& H, double e) {
  mat Hnew = H;
  ccd_update_factors(V,W,Hnew,e);
  return Hnew;
}

// Perform sequential co-ordinate descent (CCD) updates for the
// factors matrix (H) when the data matrix (W) is sparse.
//
// [[Rcpp::export]]
arma::mat ccd_update_factors_sparse_rcpp (const arma::sp_mat& V,
					  const arma::mat& W,
					  const arma::mat& H, double e) {
  mat Hnew = H;
  ccd_update_factors_sparse(V,W,Hnew,e);
  return Hnew;
}

// This does the same thing as ccd_update_factors_rcpp, except that
// Intel Threading Building Blocks (TBB) are used to update the
// loadings in parallel.
//
// [[Rcpp::depends(RcppParallel)]]
// [[Rcpp::export]]
arma::mat ccd_update_factors_parallel_rcpp (const arma::mat& V, 
					    const arma::mat& W,
					    const arma::mat& H, double e) {
  mat Hnew = H;
  ccd_factor_updater worker(V,W,Hnew,e);
  parallelFor(0,H.n_cols,worker);
  return Hnew;
}

// This does the same thing as ccd_update_factors_sparse_rcpp, except
// that Intel Threading Building Blocks (TBB) are used to update the
// loadings in parallel.
//
// [[Rcpp::export]]
arma::mat ccd_update_factors_sparse_parallel_rcpp (const arma::sp_mat& V, 
						   const arma::mat& W,
						   const arma::mat& H, 
						   double e) {
  mat Hnew = H;
  ccd_factor_updater_sparse worker(V,W,Hnew,e);
  parallelFor(0,H.n_cols,worker);
  return Hnew;
}

// Iterate the CCD updates over all columns of H, in which V is
// approximated by the matrix product W*H.
void ccd_update_factors (const mat& V, const mat& W, mat& H, double e) {
  unsigned int m = H.n_cols;
  for (unsigned int j = 0; j < m; j++)
    ccd_update_factor(V,W,H,j,e);
}

// This is the same as ccd_update_factors, except that the count data are
// stored as a sparse matrix.
void ccd_update_factors_sparse (const sp_mat& V, const mat& W, mat& H,
				double e) {
  unsigned int m = H.n_cols;
  vec sumw = sum(W,0);
  for (unsigned int j = 0; j < m; j++)
    ccd_update_factor_sparse(V,W,sumw,H,j,e);
}


<|start_filename|>docs/reference/pca_plot.html<|end_filename|>
<!-- Generated by pkgdown: do not edit by hand -->
<!DOCTYPE html>
<html lang="en">
  <head>
  <meta charset="utf-8">
<meta http-equiv="X-UA-Compatible" content="IE=edge">
<meta name="viewport" content="width=device-width, initial-scale=1.0">

<title>PCA Plot — pca_plot • fastTopics</title>


<!-- jquery -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/jquery/3.4.1/jquery.min.js" integrity="sha256-CSXorXvZcTkaix6Yvo6HppcZGetbYMGWSFlBw8HfCJo=" crossorigin="anonymous"></script>
<!-- Bootstrap -->

<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/css/bootstrap.min.css" integrity="sha256-bZLfwXAP04zRMK2BjiO8iu9pf4FbLqX6zitd+tIvLhE=" crossorigin="anonymous" />

<script src="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/js/bootstrap.min.js" integrity="sha256-nuL8/2cJ5NDSSwnKD8VqreErSWHtnEP9E7AySL+1ev4=" crossorigin="anonymous"></script>

<!-- bootstrap-toc -->
<link rel="stylesheet" href="../bootstrap-toc.css">
<script src="../bootstrap-toc.js"></script>

<!-- Font Awesome icons -->
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/all.min.css" integrity="sha256-mmgLkCYLUQbXn0B1SRqzHar6dCnv9oZFPEC1g1cwlkk=" crossorigin="anonymous" />
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/v4-shims.min.css" integrity="sha256-wZjR52fzng1pJHwx4aV2AO3yyTOXrcDW7jBpJtTwVxw=" crossorigin="anonymous" />

<!-- clipboard.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/clipboard.js/2.0.6/clipboard.min.js" integrity="sha256-inc5kl9MA1hkeYUt+EC3BhlIgyp/2jDIyBLS6k3UxPI=" crossorigin="anonymous"></script>

<!-- headroom.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/headroom.min.js" integrity="sha256-AsUX4SJE1+yuDu5+mAVzJbuYNPHj/WroHuZ8Ir/CkE0=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/jQuery.headroom.min.js" integrity="sha256-ZX/yNShbjqsohH1k95liqY9Gd8uOiE1S4vZc+9KQ1K4=" crossorigin="anonymous"></script>

<!-- pkgdown -->
<link href="../pkgdown.css" rel="stylesheet">
<script src="../pkgdown.js"></script>




<meta property="og:title" content="PCA Plot — pca_plot" />
<meta property="og:description" content="Visualize the structure of the Poisson NMF loadings or
  the multinomial topic model mixture proportions by projection onto
  principal components (PCs). plot_hexbin_plot is most useful
  for visualizing the PCs of a data set with thousands of samples, or
  more." />




<!-- mathjax -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/MathJax.js" integrity="sha256-nvJJv9wWKEm88qvoQl9ekL2J+k/RWIsaSScxxlsrv8k=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/config/TeX-AMS-MML_HTMLorMML.js" integrity="sha256-84DKXVJXs0/F8OTMzX4UR909+jtl4G7SPypPavF+GfA=" crossorigin="anonymous"></script>

<!--[if lt IE 9]>
<script src="https://oss.maxcdn.com/html5shiv/3.7.3/html5shiv.min.js"></script>
<script src="https://oss.maxcdn.com/respond/1.4.2/respond.min.js"></script>
<![endif]-->



  </head>

  <body data-spy="scroll" data-target="#toc">
    <div class="container template-reference-topic">
      <header>
      <div class="navbar navbar-default navbar-fixed-top" role="navigation">
  <div class="container">
    <div class="navbar-header">
      <button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#navbar" aria-expanded="false">
        <span class="sr-only">Toggle navigation</span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
      </button>
      <span class="navbar-brand">
        <a class="navbar-link" href="../index.html">fastTopics</a>
        <span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.5-15</span>
      </span>
    </div>

    <div id="navbar" class="navbar-collapse collapse">
      <ul class="nav navbar-nav">
        <li>
  <a href="../index.html">Home</a>
</li>
<li>
  <a href="../articles/index.html">Vignettes</a>
</li>
<li>
  <a href="../reference/index.html">Functions</a>
</li>
      </ul>
      <ul class="nav navbar-nav navbar-right">
        <li>
  <a href="https://github.com/stephenslab/fastTopics">Source</a>
</li>
      </ul>
      
    </div><!--/.nav-collapse -->
  </div><!--/.container -->
</div><!--/.navbar -->

      

      </header>

<div class="row">
  <div class="col-md-9 contents">
    <div class="page-header">
    <h1>PCA Plot</h1>
    <small class="dont-index">Source: <a href='https://github.com/stephenslab/fastTopics/blob/master/R/plots.R'><code>R/plots.R</code></a></small>
    <div class="hidden name"><code>pca_plot.Rd</code></div>
    </div>

    <div class="ref-description">
    <p>Visualize the structure of the Poisson NMF loadings or
  the multinomial topic model mixture proportions by projection onto
  principal components (PCs). <code>plot_hexbin_plot</code> is most useful
  for visualizing the PCs of a data set with thousands of samples, or
  more.</p>
    </div>

    <pre class="usage"><span class='fu'>pca_plot</span>(
  <span class='no'>fit</span>,
  <span class='no'>k</span>,
  <span class='no'>out.pca</span>,
  <span class='kw'>pcs</span> <span class='kw'>=</span> <span class='fl'>1</span>:<span class='fl'>2</span>,
  <span class='kw'>fill</span> <span class='kw'>=</span> <span class='st'>"loading"</span>,
  <span class='kw'>ggplot_call</span> <span class='kw'>=</span> <span class='no'>pca_plot_ggplot_call</span>,
  <span class='kw'>plot_grid_call</span> <span class='kw'>=</span> <span class='kw'>function</span>(<span class='no'>plots</span>) <span class='fu'><a href='https://rdrr.io/r/base/do.call.html'>do.call</a></span>(<span class='no'>plot_grid</span>, <span class='no'>plots</span>),
  <span class='no'>...</span>
)

<span class='fu'>pca_plot_ggplot_call</span>(
  <span class='no'>dat</span>,
  <span class='no'>pcs</span>,
  <span class='kw'>fill.type</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/c.html'>c</a></span>(<span class='st'>"loading"</span>, <span class='st'>"numeric"</span>, <span class='st'>"factor"</span>, <span class='st'>"none"</span>),
  <span class='no'>fill.label</span>,
  <span class='kw'>font.size</span> <span class='kw'>=</span> <span class='fl'>9</span>
)

<span class='fu'>pca_hexbin_plot</span>(
  <span class='no'>fit</span>,
  <span class='no'>out.pca</span>,
  <span class='kw'>pcs</span> <span class='kw'>=</span> <span class='fl'>1</span>:<span class='fl'>2</span>,
  <span class='kw'>bins</span> <span class='kw'>=</span> <span class='fl'>40</span>,
  <span class='kw'>breaks</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/c.html'>c</a></span>(<span class='fl'>0</span>, <span class='fl'>1</span>, <span class='fl'>10</span>, <span class='fl'>100</span>, <span class='fl'>1000</span>, <span class='fl'>Inf</span>),
  <span class='kw'>ggplot_call</span> <span class='kw'>=</span> <span class='no'>pca_hexbin_plot_ggplot_call</span>,
  <span class='no'>...</span>
)

<span class='fu'>pca_hexbin_plot_ggplot_call</span>(<span class='no'>dat</span>, <span class='no'>pcs</span>, <span class='no'>bins</span>, <span class='no'>breaks</span>, <span class='kw'>font.size</span> <span class='kw'>=</span> <span class='fl'>9</span>)</pre>

    <h2 class="hasAnchor" id="arguments"><a class="anchor" href="#arguments"></a>Arguments</h2>
    <table class="ref-arguments">
    <colgroup><col class="name" /><col class="desc" /></colgroup>
    <tr>
      <th>fit</th>
      <td><p>An object of class &#8220;poisson_nmf_fit&#8221; or
&#8220;multinom_topic_model_fit&#8221;.</p></td>
    </tr>
    <tr>
      <th>k</th>
      <td><p>The dimensions or topics selected by number or name. When
<code>fill = "loading"</code>, one plot is created per selected dimension
or topic; when <code>fill = "loading"</code> and <code>k</code> is not
specified, all dimensions or topics are plotted.</p></td>
    </tr>
    <tr>
      <th>out.pca</th>
      <td><p>A list containing the result of a principal
components analysis, typically the result of calling
<code><a href='https://rdrr.io/r/stats/prcomp.html'>prcomp</a></code> or a comparable function such as
<code>rpca</code> from the rsvd package. It should be a list object
containing, at a minimum, a list matrix or data frame
<code>out.pca$x</code> storing the rotated data. If not provided,
principal components will be computed automatically by calling
<code><a href='https://rdrr.io/r/stats/prcomp.html'>prcomp</a></code>.</p></td>
    </tr>
    <tr>
      <th>pcs</th>
      <td><p>The two principal components (PCs) to be plotted,
specified by name or number.</p></td>
    </tr>
    <tr>
      <th>fill</th>
      <td><p>The quantity to map onto the fill colour of the points
in the PCA plot. Set <code>fill = "loading"</code> to vary the fill
colour according to the loadings (or mixture proportions) of the
select topic, or topics. Alternatively, <code>fill</code> may be set to a
data vector with one entry per row of <code>fit$L</code>, in which case
these data are mapped to the fill colour of the points. When
<code>fill = "none"</code>, the fill colour is not varied.</p></td>
    </tr>
    <tr>
      <th>ggplot_call</th>
      <td><p>The function used to create the plot. Replace
<code>pca_plot_ggplot_call</code> or <code>pca_hexbin_plot_ggplot_call</code>
with your own function to customize the appearance of the plot.</p></td>
    </tr>
    <tr>
      <th>plot_grid_call</th>
      <td><p>When multiple topics (<code>k</code>) are selected,
and <code>fill = "loading"</code>, this is the function used to arrange the
plots into a grid using <code><a href='https://rdrr.io/pkg/cowplot/man/plot_grid.html'>plot_grid</a></code>. It should
be a function accepting a single argument, <code>plots</code>, a list of
<code>ggplot</code> objects.</p></td>
    </tr>
    <tr>
      <th>...</th>
      <td><p>Additional arguments passed to
<code><a href='https://rdrr.io/r/stats/prcomp.html'>prcomp</a></code>. These additional arguments are only
used if <code>out.pca</code> is not provided.</p></td>
    </tr>
    <tr>
      <th>dat</th>
      <td><p>A data frame passed as input to
<code><a href='https://ggplot2.tidyverse.org/reference/ggplot.html'>ggplot</a></code>, containing, at a minimum, the
principal components to be plotted. The data frame passed to
<code>pca_plot_ggplot_call</code> should also have a &#8220;y&#8221; column,
which is mapped to the fill colour of the points in the plot.</p></td>
    </tr>
    <tr>
      <th>fill.type</th>
      <td><p>The type of variable mapped to fill colour. The
fill colour is not varied when <code>fill.type = "none"</code>.</p></td>
    </tr>
    <tr>
      <th>fill.label</th>
      <td><p>The label used for the fill colour legend.</p></td>
    </tr>
    <tr>
      <th>font.size</th>
      <td><p>Font size used in plot.</p></td>
    </tr>
    <tr>
      <th>bins</th>
      <td><p>Number of bins used to create hexagonal 2-d
histogram. Passed as the &#8220;bins&#8221; argument to
<code><a href='https://ggplot2.tidyverse.org/reference/geom_hex.html'>stat_bin_hex</a></code>.</p></td>
    </tr>
    <tr>
      <th>breaks</th>
      <td><p>To produce the hexagonal histogram, the counts are
subdivided into intervals based on <code>breaks</code>. Passed as the
&#8220;breaks&#8221; argument to <code><a href='https://rdrr.io/r/base/cut.html'>cut</a></code>.</p></td>
    </tr>
    </table>

    <h2 class="hasAnchor" id="value"><a class="anchor" href="#value"></a>Value</h2>

    <p>A <code>ggplot</code> object.</p>
    <h2 class="hasAnchor" id="details"><a class="anchor" href="#details"></a>Details</h2>

    <p>This is a lightweight interface for rapidly producing PCA
plots from a Poisson non-negative matrix factorization or
multinomial topic model. The plotting is implemented using
ggplot2. For more control over the plot's appearance, the plot can
be customized by modifying the <code>ggplot_call</code> and
<code>plot_grid_call</code> arguments. If not provided, the principal
components are computed using <code><a href='https://rdrr.io/r/stats/prcomp.html'>prcomp</a></code>.</p>
<p>Since principal components are a <em>linear</em> projection
(&#8220;rotation&#8221;) of the data, distances between points can bee
more interpretable than nonlinear embedding methods such as t-SNE
(as visualized using <code>tsne_plot</code>) and UMAP.</p>

  </div>
  <div class="col-md-3 hidden-xs hidden-sm" id="pkgdown-sidebar">
    <nav id="toc" data-toggle="toc" class="sticky-top">
      <h2 data-toc-skip>Contents</h2>
    </nav>
  </div>
</div>


      <footer>
      <div class="copyright">
  <p>Developed by <NAME>, <NAME>, <NAME>, <NAME>.</p>
</div>

<div class="pkgdown">
  <p>Site built with <a href="https://pkgdown.r-lib.org/">pkgdown</a> 1.5.1.</p>
</div>

      </footer>
   </div>

  


  </body>
</html>




<|start_filename|>src/cost.cpp<|end_filename|>
#include "misc.h"
#include "cost.h"

using namespace arma;

// FUNCTION DEFINITIONS
// --------------------
// Compute negative log-likelihoods for assessing a topic model fit or
// quality of a non-negative matrix factorization, in which matrix X
// is approximated by matrix product A * B.
//
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::vec cost_rcpp (const arma::mat& X, const arma::mat& A,
		     const arma::mat& B, double e, bool poisson) {
  return cost(X,A,B,e,poisson);
}

// This is the same as cost_rcpp, except that X must be sparse.
//
// [[Rcpp::export]]
arma::vec cost_sparse_rcpp (const arma::sp_mat& X, const arma::mat& A,
	  		    const arma::mat& B, double e, bool poisson) {
  return cost_sparse(X,A,B,e,poisson);
}

// This is the helper function for cost_rcpp.
arma::vec cost (const mat& X, const mat& A, const mat& B, double e, 
             bool poisson) {
  unsigned int n = X.n_rows;
  unsigned int m = X.n_cols;
  vec  f(n,fill::zeros);
  vec  y(n);
  
  // Repeat for each column of X.
  for (unsigned int j = 0; j < m; j++) {

    // This is equivalent to the following R code:
    //
    //   f = f + poisson*y - X[,j]*log(y + e))
    //
    // where 
    // 
    //   y = A %*% B[,j]
    //
    y  = A * B.col(j);
    f -= X.col(j) % log(y + e);
    if (poisson)
      f += y;
  }
  
  return f;
}

// Helper function for cost_sparse_rcpp.
arma::vec cost_sparse (const sp_mat& X, const mat& A, const mat& B,
		       double e, bool poisson) {
  unsigned int n = X.n_rows;
  unsigned int m = X.n_cols;
  unsigned int i;
  vec  f(n,fill::zeros);
  vec  y(n);
  
  // Repeat for each column of X.
  for (unsigned int j = 0; j < m; j++) {

    // Initialize an iterator for the nonzero elements in the jth
    // column of X.
    sp_mat::const_col_iterator xj = X.begin_col(j);
    sp_mat::const_col_iterator xm = X.end_col(j);

    // This is equivalent to the following R code:
    //
    //   f = f + poisson*y - X[,j]*log(y + e)
    //
    // where 
    // 
    //   y = A %*% B[,j]
    //
    y = A * B.col(j);
    for(; xj != xm; ++xj) {
      i     = xj.row();
      f(i) -= (*xj) * log(y(i) + e);
    }
    if (poisson)
      f += y;
  }
  
  return f;
}


<|start_filename|>docs/reference/fit_topic_model.html<|end_filename|>
<!-- Generated by pkgdown: do not edit by hand -->
<!DOCTYPE html>
<html lang="en">
  <head>
  <meta charset="utf-8">
<meta http-equiv="X-UA-Compatible" content="IE=edge">
<meta name="viewport" content="width=device-width, initial-scale=1.0">

<title>Simple Interface for Fitting a Multinomial Topic Model — fit_topic_model • fastTopics</title>


<!-- jquery -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/jquery/3.4.1/jquery.min.js" integrity="sha256-CSXorXvZcTkaix6Yvo6HppcZGetbYMGWSFlBw8HfCJo=" crossorigin="anonymous"></script>
<!-- Bootstrap -->

<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/css/bootstrap.min.css" integrity="sha256-bZLfwXAP04zRMK2BjiO8iu9pf4FbLqX6zitd+tIvLhE=" crossorigin="anonymous" />

<script src="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/js/bootstrap.min.js" integrity="sha256-nuL8/2cJ5NDSSwnKD8VqreErSWHtnEP9E7AySL+1ev4=" crossorigin="anonymous"></script>

<!-- bootstrap-toc -->
<link rel="stylesheet" href="../bootstrap-toc.css">
<script src="../bootstrap-toc.js"></script>

<!-- Font Awesome icons -->
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/all.min.css" integrity="sha256-mmgLkCYLUQbXn0B1SRqzHar6dCnv9oZFPEC1g1cwlkk=" crossorigin="anonymous" />
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/v4-shims.min.css" integrity="sha256-wZjR52fzng1pJHwx4aV2AO3yyTOXrcDW7jBpJtTwVxw=" crossorigin="anonymous" />

<!-- clipboard.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/clipboard.js/2.0.6/clipboard.min.js" integrity="sha256-inc5kl9MA1hkeYUt+EC3BhlIgyp/2jDIyBLS6k3UxPI=" crossorigin="anonymous"></script>

<!-- headroom.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/headroom.min.js" integrity="sha256-AsUX4SJE1+yuDu5+mAVzJbuYNPHj/WroHuZ8Ir/CkE0=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/jQuery.headroom.min.js" integrity="sha256-ZX/yNShbjqsohH1k95liqY9Gd8uOiE1S4vZc+9KQ1K4=" crossorigin="anonymous"></script>

<!-- pkgdown -->
<link href="../pkgdown.css" rel="stylesheet">
<script src="../pkgdown.js"></script>




<meta property="og:title" content="Simple Interface for Fitting a Multinomial Topic Model — fit_topic_model" />
<meta property="og:description" content="Fits a multinomial topic model to the count data,
  hiding most of the complexities of model fitting. The default
  optimization settings used here are intended to work well in a wide
  range of data sets, although some fine-tuning may be needed for
  more difficult cases. For full control, use fit_poisson_nmf." />




<!-- mathjax -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/MathJax.js" integrity="sha256-nvJJv9wWKEm88qvoQl9ekL2J+k/RWIsaSScxxlsrv8k=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/config/TeX-AMS-MML_HTMLorMML.js" integrity="sha256-84DKXVJXs0/F8OTMzX4UR909+jtl4G7SPypPavF+GfA=" crossorigin="anonymous"></script>

<!--[if lt IE 9]>
<script src="https://oss.maxcdn.com/html5shiv/3.7.3/html5shiv.min.js"></script>
<script src="https://oss.maxcdn.com/respond/1.4.2/respond.min.js"></script>
<![endif]-->



  </head>

  <body data-spy="scroll" data-target="#toc">
    <div class="container template-reference-topic">
      <header>
      <div class="navbar navbar-default navbar-fixed-top" role="navigation">
  <div class="container">
    <div class="navbar-header">
      <button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#navbar" aria-expanded="false">
        <span class="sr-only">Toggle navigation</span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
      </button>
      <span class="navbar-brand">
        <a class="navbar-link" href="../index.html">fastTopics</a>
        <span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.5-52</span>
      </span>
    </div>

    <div id="navbar" class="navbar-collapse collapse">
      <ul class="nav navbar-nav">
        <li>
  <a href="../index.html">Home</a>
</li>
<li>
  <a href="../articles/index.html">Vignettes</a>
</li>
<li>
  <a href="../reference/index.html">Functions</a>
</li>
      </ul>
      <ul class="nav navbar-nav navbar-right">
        <li>
  <a href="https://github.com/stephenslab/fastTopics">Source</a>
</li>
      </ul>
      
    </div><!--/.nav-collapse -->
  </div><!--/.container -->
</div><!--/.navbar -->

      

      </header>

<div class="row">
  <div class="col-md-9 contents">
    <div class="page-header">
    <h1>Simple Interface for Fitting a Multinomial Topic Model</h1>
    <small class="dont-index">Source: <a href='https://github.com/stephenslab/fastTopics/blob/master/R/fit_topic_model.R'><code>R/fit_topic_model.R</code></a></small>
    <div class="hidden name"><code>fit_topic_model.Rd</code></div>
    </div>

    <div class="ref-description">
    <p>Fits a multinomial topic model to the count data,
  hiding most of the complexities of model fitting. The default
  optimization settings used here are intended to work well in a wide
  range of data sets, although some fine-tuning may be needed for
  more difficult cases. For full control, use <code>fit_poisson_nmf</code>.</p>
    </div>

    <pre class="usage"><span class='fu'>fit_topic_model</span>(
  <span class='no'>X</span>,
  <span class='no'>k</span>,
  <span class='kw'>numiter.main</span> <span class='kw'>=</span> <span class='fl'>100</span>,
  <span class='kw'>numiter.refine</span> <span class='kw'>=</span> <span class='fl'>100</span>,
  <span class='kw'>method.main</span> <span class='kw'>=</span> <span class='st'>"em"</span>,
  <span class='kw'>method.refine</span> <span class='kw'>=</span> <span class='st'>"scd"</span>,
  <span class='kw'>init.method</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/c.html'>c</a></span>(<span class='st'>"topicscore"</span>, <span class='st'>"random"</span>),
  <span class='kw'>control.init</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/list.html'>list</a></span>(),
  <span class='kw'>control.main</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/list.html'>list</a></span>(<span class='kw'>numiter</span> <span class='kw'>=</span> <span class='fl'>4</span>),
  <span class='kw'>control.refine</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/list.html'>list</a></span>(<span class='kw'>numiter</span> <span class='kw'>=</span> <span class='fl'>4</span>, <span class='kw'>extrapolate</span> <span class='kw'>=</span> <span class='fl'>TRUE</span>),
  <span class='kw'>verbose</span> <span class='kw'>=</span> <span class='fu'><a href='https://rdrr.io/r/base/c.html'>c</a></span>(<span class='st'>"progressbar"</span>, <span class='st'>"detailed"</span>, <span class='st'>"none"</span>)
)</pre>

    <h2 class="hasAnchor" id="arguments"><a class="anchor" href="#arguments"></a>Arguments</h2>
    <table class="ref-arguments">
    <colgroup><col class="name" /><col class="desc" /></colgroup>
    <tr>
      <th>X</th>
      <td><p>The n x m matrix of counts; all entries of X should be
non-negative. It can be a sparse matrix (class <code>"dgCMatrix"</code>)
or dense matrix (class <code>"matrix"</code>).</p></td>
    </tr>
    <tr>
      <th>k</th>
      <td><p>The number of topics. Must be 2 or greater.</p></td>
    </tr>
    <tr>
      <th>numiter.main</th>
      <td><p>Number of updates of the factors and loadings
to perform in the main model fitting step. Should be 1 or more.</p></td>
    </tr>
    <tr>
      <th>numiter.refine</th>
      <td><p>Number of updates of the factors and loadings
to perform in the model refinement step.</p></td>
    </tr>
    <tr>
      <th>method.main</th>
      <td><p>The method to use for updating the factors and
loadings in the main model fitting step. Passed as argument
"method" to <code><a href='fit_poisson_nmf.html'>fit_poisson_nmf</a></code>.</p></td>
    </tr>
    <tr>
      <th>method.refine</th>
      <td><p>The method to use for updating the factors in
evthe model refinement step. Passed as argument "method"
to <code><a href='fit_poisson_nmf.html'>fit_poisson_nmf</a></code>.</p></td>
    </tr>
    <tr>
      <th>init.method</th>
      <td><p>The method used to initialize the factors and
loadings. See <code><a href='fit_poisson_nmf.html'>init_poisson_nmf</a></code> for details.</p></td>
    </tr>
    <tr>
      <th>control.init</th>
      <td><p>A list of parameters controlling the behaviour
of the optimization and Topic SCORE method in the call to
<code>init_poisson_nmf</code>. This is passed as argument "control" to
<code>init_poisson_nmf</code>.</p></td>
    </tr>
    <tr>
      <th>control.main</th>
      <td><p>A list of parameters controlling the behaviour
of the optimization in the main model fitting step. This is passed
as argument "control" to <code>fit_poisson_nmf</code>.</p></td>
    </tr>
    <tr>
      <th>control.refine</th>
      <td><p>A list of parameters controlling the
behaviour of the of the optimization algorithm in the model
refinement step. This is passed as argument "control" to
<code>fit_poisson_nmf</code>.</p></td>
    </tr>
    <tr>
      <th>verbose</th>
      <td><p>When <code>verbose = "progressbar"</code> or <code>verbose
= "detailed"</code>, information about the progress of the model fitting
is printed to the console. See <code><a href='fit_poisson_nmf.html'>fit_poisson_nmf</a></code> for
more information.</p></td>
    </tr>
    </table>

    <h2 class="hasAnchor" id="value"><a class="anchor" href="#value"></a>Value</h2>

    <p>A multinomial topic model fit; see
  <code><a href='poisson2multinom.html'>poisson2multinom</a></code> and <code><a href='fit_poisson_nmf.html'>fit_poisson_nmf</a></code>
  for details. Note that outputted likelihoods and deviances in
  <code>progress</code> are evaluated with respect to the equivalent
  Poisson NMF model.</p>
    <h2 class="hasAnchor" id="details"><a class="anchor" href="#details"></a>Details</h2>

    <p>With the default settings, the model fitting is
accomplished in four steps: (1) initialize the Poisson NMF model
fit (<code>init_poisson_nmf</code>); (2) perform the main model fitting
step by running 100 EM updates using <code>fit_poisson_nmf</code>; (3)
refine the fit by running 100 extrapolated SCD updates, again using
<code>fit_poisson_nmf</code>; and (4) recover the multinomial topic model
by calling <code>poisson2multinom</code>.</p>
<p>This two-stage fitting approach is based on our findings that the
EM algorithm initially makes rapid progress toward a solution, but
its convergence slows considerably as the iterates approach a
solution. Close to a solution, we have found that other algorithms
make much more rapid progress than EM; in particularly, we founr
that the extrapolated SCD updates usually performed best). For
larger data sets, more updates in the main model fitting and
refinement steps may be needed to obtain a good fit.</p>
<p>For larger data sets, more than 200 updates may be needed to obtain
a good fit.</p>
    <h2 class="hasAnchor" id="references"><a class="anchor" href="#references"></a>References</h2>

    <p><NAME>., <NAME>. and <NAME>. (2017). Visualizing the
structure of RNA-seq expression data using grade of membership
models. <em>PLoS Genetics</em> <b>13</b>, e1006599.</p>
<p><NAME>., <NAME>. and <NAME>. (2003). Latent Dirichlet
allocation. <em>Journal of Machine Learning Research</em> <b>3</b>,
993-1022.</p>
<p><NAME>. (1999). Probabilistic latent semantic indexing. In
<em>Proceedings of the 22nd International ACM SIGIR Conference</em>,
50-57. doi:10.1145/312624.312649</p>
    <h2 class="hasAnchor" id="see-also"><a class="anchor" href="#see-also"></a>See also</h2>

    <div class='dont-index'><p><code><a href='fit_poisson_nmf.html'>init_poisson_nmf</a></code>,
  <code><a href='fit_poisson_nmf.html'>fit_poisson_nmf</a></code>,
  <code><a href='poisson2multinom.html'>poisson2multinom</a></code>,
  <code><a href='fit_multinom_model.html'>fit_multinom_model</a></code></p></div>

    <h2 class="hasAnchor" id="examples"><a class="anchor" href="#examples"></a>Examples</h2>
    <pre class="examples"><div class='input'><span class='fu'><a href='https://rdrr.io/r/base/library.html'>library</a></span>(<span class='no'>Matrix</span>)
<span class='fu'><a href='https://rdrr.io/r/base/Random.html'>set.seed</a></span>(<span class='fl'>1</span>)
<span class='no'>X</span> <span class='kw'>&lt;-</span> <span class='fu'><a href='simulate_count_data.html'>simulate_count_data</a></span>(<span class='fl'>80</span>,<span class='fl'>100</span>,<span class='kw'>k</span> <span class='kw'>=</span> <span class='fl'>3</span>,<span class='kw'>sparse</span> <span class='kw'>=</span> <span class='fl'>TRUE</span>)$<span class='no'>X</span>
<span class='no'>fit</span> <span class='kw'>&lt;-</span> <span class='fu'>fit_topic_model</span>(<span class='no'>X</span>,<span class='kw'>k</span> <span class='kw'>=</span> <span class='fl'>3</span>)</div><div class='output co'>#&gt; Initializing factors using Topic SCORE algorithm.
#&gt; Initializing loadings by running 10 SCD updates.</div><div class='output co'>#&gt; <span class='message'>Using 4 RcppParallel threads.</span></div><div class='output co'>#&gt; <span class='message'>Using 4 RcppParallel threads.</span></div><div class='output co'>#&gt; Fitting rank-3 Poisson NMF to 80 x 100 sparse matrix.
#&gt; Running 100 EM updates, without extrapolation (fastTopics 0.5-52).
#&gt; Refining model fit.</div><div class='output co'>#&gt; <span class='message'>Using 4 RcppParallel threads.</span></div><div class='output co'>#&gt; Fitting rank-3 Poisson NMF to 80 x 100 sparse matrix.
#&gt; Running 100 SCD updates, with extrapolation (fastTopics 0.5-52).</div><div class='input'><span class='fu'><a href='https://rdrr.io/r/base/print.html'>print</a></span>(<span class='fu'><a href='https://rdrr.io/r/base/summary.html'>summary</a></span>(<span class='no'>fit</span>))</div><div class='output co'>#&gt; Model overview:
#&gt;   Number of data rows, n: 80
#&gt;   Number of data cols, m: 100
#&gt;   Rank/number of topics, k: 3
#&gt; Evaluation of model fit (200 updates performed):
#&gt;   Poisson NMF log-likelihood: -8.311490191185e+03
#&gt;   Multinomial topic model log-likelihood: -8.068320081832e+03
#&gt;   Poisson NMF deviance: +7.710268688898e+03
#&gt;   Max KKT residual: +9.610230e-08</div><div class='output co'>#&gt; <span class='message'>Set show.size.factors = TRUE, show.mixprops = TRUE and/or show.topic.reps = TRUE in print(...) to show more informataion</span></div><div class='input'>
</div></pre>
  </div>
  <div class="col-md-3 hidden-xs hidden-sm" id="pkgdown-sidebar">
    <nav id="toc" data-toggle="toc" class="sticky-top">
      <h2 data-toc-skip>Contents</h2>
    </nav>
  </div>
</div>


      <footer>
      <div class="copyright">
  <p>Developed by <NAME>, <NAME>, <NAME>, <NAME>.</p>
</div>

<div class="pkgdown">
  <p>Site built with <a href="https://pkgdown.r-lib.org/">pkgdown</a> 1.5.1.</p>
</div>

      </footer>
   </div>

  


  </body>
</html>




<|start_filename|>src/pnmfem.cpp<|end_filename|>
#include <RcppParallel.h>
#include "misc.h"
#include "poismix.h"

using namespace arma;

// INLINE FUNCTION DEFINITIONS
// ---------------------------
// Perform one or several EM updates for a single column of the k x m
// factors matrix, F, when X is a dense matrix.
inline vec pnmfem_update_factor (const mat& X, const mat& F, const mat& L1,
				 const vec& u, mat& P, unsigned int j, 
				 unsigned int numiter) {
  vec f = F.col(j);
  poismixem(L1,u,X.col(j),f,P,numiter);
  return f;
}

// Perform one or several EM updates for a single column of the k x m
// factors matrix, F, when X is a sparse matrix.
inline vec pnmfem_update_factor_sparse (const sp_mat& X, const mat& F,
					const mat& L1, const vec& u, 
					unsigned int j, unsigned int numiter) {
  vec          x = nonzeros(X.col(j));
  unsigned int n = x.n_elem;
  vec          f = F.col(j);
  uvec         i(n);
  getcolnonzeros(X,i,j);
  poismixem(L1,u,x,i,f,numiter);
  return f;
}

// CLASS DEFINITIONS
// -----------------
// This class is used to implement multithreaded computation of the
// factor updates in pnmfem_update_factors_parallel_rcpp.
struct pnmfem_factor_updater : public RcppParallel::Worker {
  const mat&   X;
  const mat&   F;
  mat          L1;
  vec          u;
  mat&         Fnew;
  const vec&   j;
  unsigned int numiter;

  // This is used to create a pnmfem_factor_updater object.
  pnmfem_factor_updater (const mat& X, const mat& F, const mat& L, 
			 mat& Fnew, const vec& j, unsigned int numiter) :
    X(X), F(F), L1(L), u(L.n_cols), Fnew(Fnew), j(j), numiter(numiter) { 
    u = sum(L,0);
    normalizecols(L1);
  };

  // This function updates the factors for a given range of columns.
  void operator() (std::size_t begin, std::size_t end) {
    mat P = L1;
    for (unsigned int i = begin; i < end; i++)
      Fnew.col(j(i)) = pnmfem_update_factor(X,F,L1,u,P,j(i),numiter);
  }
};

// This class is used to implement multithreaded computation of the
// factor updates in pnmfem_update_factors_sparse_parallel_rcpp.
struct pnmfem_factor_updater_sparse : public RcppParallel::Worker {
  const sp_mat& X;
  const mat&    F;
  mat           L1;
  vec           u;
  mat&          Fnew;
  const vec&    j;
  unsigned int  numiter;

  // This is used to create a pnmfem_factor_updater object.
  pnmfem_factor_updater_sparse (const sp_mat& X, const mat& F, const mat& L, 
				mat& Fnew, const vec& j, 
				unsigned int numiter) :
    X(X), F(F), L1(L), u(L.n_cols), Fnew(Fnew), j(j), numiter(numiter) {
    u = sum(L,0);
    normalizecols(L1);
  };

  // This function updates the factors for a given range of columns.
  void operator() (std::size_t begin, std::size_t end) {
    for (unsigned int i = begin; i < end; i++)
      Fnew.col(j(i)) = pnmfem_update_factor_sparse(X,F,L1,u,j(i),numiter);
  }
};

// FUNCTION DEFINITIONS
// --------------------
// Perform one or several EM updates for the factors matrix, F, in
// which the matrix X is approximated by L*F. Input "k" specifies
// which columns of F are updated. Input "numiter" specifies the
// number of EM updates to perform.
//
// Note that, unlike most other functions implemented in this package,
// the factors matrix F should be a k x m matrix, where k is the
// number factors, or "topics". This is done for ease of
// implementation, and for speed, because the computation is performed
// on the m columns of F, and the matrix is stored columnwise.
//
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat pnmfem_update_factors_rcpp (const arma::mat& X, const arma::mat& F,
				      const arma::mat& L, const arma::vec& j,
				      double numiter) {
  unsigned int n = j.n_elem;
  vec  u    = sum(L,0);
  mat  L1   = L;
  mat  P    = L;
  mat  Fnew = F;
  normalizecols(L1);
  for (unsigned int i = 0; i < n; i++)
    Fnew.col(j(i)) = pnmfem_update_factor(X,F,L1,u,P,j(i),numiter);
  return Fnew;
}

// This does the same thing as pnmfem_update_factors_rcpp, except that
// X is a sparse matrix. See pnmfem_update_factors_rcpp for details.
//
// [[Rcpp::export]]
arma::mat pnmfem_update_factors_sparse_rcpp (const arma::sp_mat& X,
					     const arma::mat& F,
					     const arma::mat& L,
					     const arma::vec& j,
					     double numiter) {
  unsigned int n = j.n_elem;
  vec  u    = sum(L,0);
  mat  L1   = L;
  mat  Fnew = F;
  normalizecols(L1);
  for (unsigned int i = 0; i < n; i++)
    Fnew.col(j(i)) = pnmfem_update_factor_sparse(X,F,L1,u,j(i),numiter);
  return Fnew;
}

// This does the same thing as pnmfem_update_factors_rcpp, except that
// Intel Threading Building Blocks (TBB) are used to update the
// factors in parallel.
//
// [[Rcpp::depends(RcppParallel)]]
// [[Rcpp::export]]
arma::mat pnmfem_update_factors_parallel_rcpp (const arma::mat& X,
					       const arma::mat& F,
					       const arma::mat& L,
					       const arma::vec& j,
					       double numiter) {
  mat Fnew = F;
  pnmfem_factor_updater worker(X,F,L,Fnew,j,numiter);
  parallelFor(0,j.n_elem,worker);
  return Fnew;
}

// This does the same thing as pnmfem_update_factors_sparse_rcpp,
// except that Intel Threading Building Blocks (TBB) are used to
// update the factors in parallel.
//
// [[Rcpp::export]]
arma::mat pnmfem_update_factors_sparse_parallel_rcpp (const arma::sp_mat& X,
						      const arma::mat& F,
						      const arma::mat& L,
						      const arma::vec& j,
						      double numiter) {
  mat Fnew = F;
  pnmfem_factor_updater_sparse worker(X,F,L,Fnew,j,numiter);
  parallelFor(0,j.n_elem,worker);
  return Fnew;
}


<|start_filename|>src/poismix.cpp<|end_filename|>
#include "poismix.h"
#include "misc.h"
#include "mixem.h"

using namespace arma;

// INLINE FUNCTION DEFINITIONS
// ---------------------------
// This is the same as compute_grad_scd (see below), except that u =
// sum(L[,j]), w must contain only the nonzero values, and vectors l
// and Lx should contain only the entries associated with these
// nonzero values (vectors l, w and Lx should be the same length). 
// This is used by scd_kl_update.
inline void compute_grad_scd_sparse (const vec& l, const vec& w, const vec& Lx,
				     double u, double x, double& g, double& h,
				     vec& r, double e) {
  r = l/(Lx + e);
  h = dot(square(r),w);
  g = dot(r,w) - u + h*x;
}

// Compute the first-order and second-order partial derivatives of the
// loss function (i.e., the negative Poisson log-likelihood) with
// respect to the jth mixture weight, x[j], where: input argument l =
// L[,j] is the jth column of the n x m data matrix, L; input w is a
// vector of length n containing a count or "pseudocount" associated
// with each row of L; x is the initial estimate of the mixture
// weights; and Lx = L*x. Input argument u is a vector storing an
// intermediate result of the same size as w. The outputs are g = h*x
// - d and h, where d and h are the first-order and second-order
// partial derivatives of the loss function. This is used by
// scd_kl_update.
inline void compute_grad_scd (const vec& l, const vec& w, const vec& Lx,
			      double x, double& g, double& h, vec& r,
			      double e) {
  compute_grad_scd_sparse(l,w,Lx,sum(l),x,g,h,r,e);
}

// Given g = h*x - d and h, where d and h are the first-order and
// second-order partial derivatives of the loss function (negative
// Poisson log-likelihood), compute the SCD update of the mixture
// weight, x[j]. This is used by scd_kl_update.
inline double project_iterate_scd (double g, double h, double e) {
  double x = g/(h + e);
  x = maximum(x,0);
  return x;
}

// Compute the first-order (g) and second-order (h) partial
// derivatives of the loss function (i.e., the negative Poisson
// log-likelihood) with respect to the jth mixture weight, x[j]. See
// compute_grad_scd for an explanation of the inputs. This is used by
// ccd_kl_update.
inline void compute_grad_ccd (const vec& l, const vec& w, const vec& Lx,
			      double x, double& g, double& h, double e) {
  unsigned int n = l.n_elem;
  double       t;
  g = 0;
  h = 0;
  for (unsigned int i = 0; i < n; i++) {
    t  = w(i)/(Lx(i) + e);
    g += l(i)*(1 - t);
    h += l(i)*l(i)*t/(Lx(i) + e);
  }
}

// This is the same as compute_grad_ccd, except that u = sum(L[,j]), w
// must contain only the nonzero values, and vectors l and Lx should
// contain only the entries associated with these nonzero values
// (vectors l, w and Lx should be the same length). This is used by
// ccd_kl_update.
inline void compute_grad_ccd_sparse (const vec& l, const vec& w, const vec& Lx,
				     double u, double x, double& g, double& h,
				     double e) {
  unsigned int n = w.n_elem;
  double       t;
  g = u;
  h = 0;
  for (unsigned int i = 0; i < n; i++) {
    t  = w(i)/(Lx(i) + e);
    g -= l(i)*t;
    h += l(i)*l(i)*t/(Lx(i) + e);
  }
}

// Given the first-order (g) and second-order (h) partial derivatives
// of the loss function (negative Poisson log-likelihood), compute the
// CCD update of the mixture weight, w[j]. This is used by ccd_kl_update.
inline double project_iterate_ccd (double x, double g, double h, double e) {
  double y;
  y = x - g/h + e;
  y = maximum(y,e);
  return y;
}

// FUNCTION DEFINITIONS
// --------------------
// This is mainly used to test the first variant of the poismixem C++
// function.
//
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::vec poismixem_rcpp (const arma::mat& L, const arma::vec& w,
			  const arma::vec& x0, unsigned int numiter) {
  return poismixem(L,w,x0,numiter);
}

// This is mainly used to test the second variant of the poismixem C++
// interface.
//
// [[Rcpp::export]]
arma::vec poismixem2_rcpp (const arma::mat& L1, const arma::vec& w,
			   const arma::vec& u, const arma::vec& x0,
			   unsigned int numiter) {
  vec x = x0;
  mat P = L1;
  poismixem(L1,u,w,x,P,numiter);
  return x;
}

// This is mainly used to test the third variant of the poismixem C++
// interface.
//
// [[Rcpp::export]]
arma::vec poismixem3_rcpp (const arma::mat& L1, const arma::vec& w,
			   const arma::vec& u, const arma::uvec& i,
			   const arma::vec& x0, unsigned int numiter) {
  vec x = x0;
  poismixem(L1,u,w,i,x,numiter);
  return x;
}

// This is mainly used to test the first variant of the C++ function
// scd_kl_update.
//
// [[Rcpp::export]]
arma::vec scd_kl_update_rcpp (const arma::mat& L, const arma::vec& w, 
			      const arma::vec& x0, unsigned int numiter, double e) {
  return scd_kl_update(L,w,x0,numiter,e);
}

// This is mainly used to test the second variant of the C++ function
// scd_kl_update.
//
// [[Rcpp::export]]
arma::vec scd_kl_update2_rcpp (const arma::mat& L, const arma::vec& u,
			       const arma::vec& w, const arma::vec& x0,
			       unsigned int numiter, double e) {
  return scd_kl_update(L,u,w,x0,numiter,e);
}

// This is mainly used to test the first variant of the C++ function
// ccd_kl_update.
//
// [[Rcpp::export]]
arma::vec ccd_kl_update_rcpp (const arma::mat& L, const arma::vec& w,
			      const arma::vec& x0, unsigned int numiter, 
			      double e) {
  vec x = x0;
  for (int iter = 0; iter < numiter; iter++)
    x = ccd_kl_update(L,w,x,e);
  return x;
}

// This is mainly used to test the second variant of the C++ function
// ccd_kl_update.
//
// [[Rcpp::export]]
arma::vec ccd_kl_update2_rcpp (const arma::mat& L, const arma::vec& u,
			       const arma::vec& w, const arma::vec& x0,
			       unsigned int numiter, double e) {
  vec x = x0;
  for (int iter = 0; iter < numiter; iter++) 
    x = ccd_kl_update(L,u,w,x,e);
  return x;
}

// Compute a maximum-likelihood estimate (MLE) of the mixture weights
// in a Poisson mixture model by iterating the multinomial mixture
// model EM updates for a fixed number of iterations.
//
// Input argument L is an n x m matrix with non-negative entries;
// input w is a vector of length n containing a count or "pseudocount"
// associated with each row of L; input argument x0 is the initial
// estimate of the mixture weights; and input "numiter" specifies the
// number of EM updates to perform.
//
// The return value is a vector of length m containing the updated
// mixture weights.
vec poismixem (const mat& L, const vec& w, const vec& x0, 
	       unsigned int numiter) {
  mat L1 = L;
  mat P  = L;
  vec u  = sum(L,0);
  vec x  = x0;
  normalizecols(L1);
  poismixem(L1,u,w,x,P,numiter);
  return x;
}

// Use this variant of poismixem if you plan on calling poismixem
// multiple times with the same matrix, L. In this call, input u should
// contain the column sums, u = colSums(L), in which L is the matrix
// prior to normalization, and matrix L1 is the normalized version of
// L in which each column sums to 1; that is, L1 = normalize.cols(L).
// Input P should be a matrix of the same size as L1.
//
// Note that in this variant of poismixem, L1 and w do not need to
// contain all the data; once u is correctly calculated, any rows of
// L1 associated with zero weights have no effect, so only the vector
// of nonzero weights w, and the rows of L1 associated with those
// weights, need to be supplied.
void poismixem (const mat& L1, const vec& u, const vec& w, vec& x, mat& P, 
		unsigned int numiter) {
  double s = sum(w);
  
  // Recover the mixture proportions of the multinomial mixture model
  // from the mixture weights of the Poisson mixture model. 
  x %= u;
  x /= sum(x);

  // Perform one or more EM updates for the multinomial mixture model.
  mixem(L1,w,x,P,numiter);

  // Recover the mixture weights of the Poisson mixture model from the
  // mixture weights of the multinomial mixture model.
  x *= s;
  x /= u;
}

// This third poismixem interface is similar to the second, except
// that the indices of the nonzero counts are supplied in vector i,
// and w must contain only the nonzero counts (vectors i and w should
// be the same length).
//
// To illustrate, in this example x1 and x2 should be the same after
// running this R code:
//
//   i  <- which(w > 0)
//   x1 <- poismixem2_rcpp(L1,w,u,x0,numiter)
//   x2 <- poismixem3_rcpp(L1,w[i],u,i-1,x0,numiter)
//
void poismixem (const mat& L1, const vec& u, const vec& w, const uvec& i,
		vec& x, unsigned int numiter) {
  unsigned int n = i.n_elem;
  unsigned int m = x.n_elem;
  mat  P(n,m);
  poismixem(L1.rows(i),u,w,x,P,numiter);
}

// Compute a maximum-likelihood estimate (MLE) of the mixture weights
// in a Poisson mixture model by iterating coordinate-wise updates for
// a fixed number of iterations.
//
// Input argument L is an n x m matrix with non-negative entries;
// input w is a vector of length n containing a count or "pseudocount"
// associated with each row of L; input argument x0 is the initial
// estimate of the mixture weights; and input "numiter" specifies the
// number of updates (full passes) to perform.
//
// The return value is a vector of length m containing the updated
// mixture weights.
//  
// This function implements the core part of the sequential
// co-ordinate descent (SCD) algorithm. This C++ code is adapted from
// the code by <NAME> and <NAME>, which is available for
// download at https://github.com/linxihui/NNLM.
vec scd_kl_update (const mat& L, const vec& w, const vec& x0,
		   unsigned int numiter, double e) {
  unsigned int n = L.n_rows;
  unsigned int m = L.n_cols;
  vec    x  = x0;
  vec    Lx = L * x;
  vec    l(n);
  vec    r(n);
  double h, g, xj, xjnew;
  for (unsigned int iter = 0; iter < numiter; iter++)
    for (unsigned int j = 0; j < m; j++) {
      l     = L.col(j);
      xj    = x(j);
      compute_grad_scd(l,w,Lx,xj,g,h,r,e);
      xjnew = project_iterate_scd(g,h,e);
      Lx   += (xjnew - xj) * l;
      x(j)  = xjnew;
    }
  return x;
}

// In this alternative implementation of the SCD updates, w must
// contain only the nonzero counts, and L must contain only the rows
// of the matrix associated with the nonzero counts (vectors i and w
// should be the same length, and the number of rows in L should be
// the same as the lengths of vectors i and w).
vec scd_kl_update (const mat& L, const vec& u, const vec& w,
		   const vec& x0, unsigned int numiter, double e) {
  unsigned int n = w.n_elem;
  unsigned int m = L.n_cols;
  vec    x  = x0;
  vec    Lx = L * x;
  vec    l(n);
  vec    r(n);
  double h, g, xj, xjnew;
  for (unsigned int iter = 0; iter < numiter; iter++)
    for (unsigned int j = 0; j < m; j++) {
      l     = L.col(j);
      xj    = x(j);
      compute_grad_scd_sparse(l,w,Lx,u(j),xj,g,h,r,e);
      xjnew = project_iterate_scd(g,h,e);
      Lx   += (xjnew - xj) * l;
      x(j)  = xjnew;
    }
  return x;
}

// Implements a single iteration of the cyclic co-ordinate descent (CCD)
// algorithm. See scd_kl_update for an explanation of the inputs.
//
// This implementation is adapted from the C++ code developed by
// <NAME> and <NAME>, which is available for download
// at www.cs.utexas.edu/~cjhsieh/nmf.
vec ccd_kl_update (const mat& L, const vec& w, const vec& x0, double e) {
  unsigned int n = L.n_rows;
  unsigned int m = L.n_cols;
  vec    x  = x0;
  vec    Lx = L * x;
  vec    l(n);
  double g, h, xj, xjnew;
  for (unsigned int j = 0; j < m; j++) {
    l     = L.col(j);
    xj    = x(j);
    compute_grad_ccd(l,w,Lx,xj,g,h,e);
    xjnew = project_iterate_ccd(xj,g,h,e);
    Lx   += (xjnew - xj) * l;
    x(j)  = xjnew;
  }

  return x;
}

// In this alternative implementation of the CCD updates, w must
// contain only the nonzero counts, and L must contain only the rows
// of the matrix associated with the nonzero counts (vectors i and w
// should be the same length, and the number of rows in L should be
// the same as the lengths of vectors i and w).
vec ccd_kl_update (const mat& L, const vec& u, const vec& w,
		   const vec& x0, double e) {
  unsigned int n = w.n_elem;
  unsigned int m = L.n_cols;
  vec    x  = x0;
  vec    Lx = L * x;
  vec    l(n);
  double h, g, xj, xjnew;
  for (unsigned int j = 0; j < m; j++) {
    l     = L.col(j);
    xj    = x(j);
    compute_grad_ccd_sparse(l,w,Lx,u(j),xj,g,h,e);
    xjnew = project_iterate_ccd(xj,g,h,e);
    Lx   += (xjnew - xj) * l;
    x(j)  = xjnew;
  }
  return x;
}


<|start_filename|>src/mixem.h<|end_filename|>
#ifndef INCLUDE_MIXEM
#define INCLUDE_MIXEM

#include <RcppArmadillo.h>

// FUNCTION DECLARATIONS
// ---------------------
arma::vec mixem (const arma::mat& L, const arma::vec& w, const arma::vec& x0,
		 unsigned int numiter);

void mixem (const arma::mat& L1, const arma::vec& w, arma::vec& x,
	    arma::mat& P, unsigned int numiter);

void mixture_one_nonzero (const arma::mat& L1, unsigned int i, arma::vec& x);

#endif


<|start_filename|>docs/reference/tsne_from_topics.html<|end_filename|>
<!-- Generated by pkgdown: do not edit by hand -->
<!DOCTYPE html>
<html lang="en">
  <head>
  <meta charset="utf-8">
<meta http-equiv="X-UA-Compatible" content="IE=edge">
<meta name="viewport" content="width=device-width, initial-scale=1.0">

<title>t-SNE from Poisson NMF or Multinomial Topic Model — tsne_from_topics • fastTopics</title>


<!-- jquery -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/jquery/3.4.1/jquery.min.js" integrity="sha256-CSXorXvZcTkaix6Yvo6HppcZGetbYMGWSFlBw8HfCJo=" crossorigin="anonymous"></script>
<!-- Bootstrap -->

<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/css/bootstrap.min.css" integrity="sha256-bZLfwXAP04zRMK2BjiO8iu9pf4FbLqX6zitd+tIvLhE=" crossorigin="anonymous" />

<script src="https://cdnjs.cloudflare.com/ajax/libs/twitter-bootstrap/3.4.1/js/bootstrap.min.js" integrity="sha256-nuL8/2cJ5NDSSwnKD8VqreErSWHtnEP9E7AySL+1ev4=" crossorigin="anonymous"></script>

<!-- bootstrap-toc -->
<link rel="stylesheet" href="../bootstrap-toc.css">
<script src="../bootstrap-toc.js"></script>

<!-- Font Awesome icons -->
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/all.min.css" integrity="sha256-mmgLkCYLUQbXn0B1SRqzHar6dCnv9oZFPEC1g1cwlkk=" crossorigin="anonymous" />
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.12.1/css/v4-shims.min.css" integrity="sha256-wZjR52fzng1pJHwx4aV2AO3yyTOXrcDW7jBpJtTwVxw=" crossorigin="anonymous" />

<!-- clipboard.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/clipboard.js/2.0.6/clipboard.min.js" integrity="sha256-inc5kl9MA1hkeYUt+EC3BhlIgyp/2jDIyBLS6k3UxPI=" crossorigin="anonymous"></script>

<!-- headroom.js -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/headroom.min.js" integrity="sha256-AsUX4SJE1+yuDu5+mAVzJbuYNPHj/WroHuZ8Ir/CkE0=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/headroom/0.11.0/jQuery.headroom.min.js" integrity="sha256-ZX/yNShbjqsohH1k95liqY9Gd8uOiE1S4vZc+9KQ1K4=" crossorigin="anonymous"></script>

<!-- pkgdown -->
<link href="../pkgdown.css" rel="stylesheet">
<script src="../pkgdown.js"></script>




<meta property="og:title" content="t-SNE from Poisson NMF or Multinomial Topic Model — tsne_from_topics" />
<meta property="og:description" content="Computes a low-dimensional nonlinear embededding of
  the data from the estimated loadings or mixture proportions using
  the t-SNE nonlinear dimensionality reduction method." />




<!-- mathjax -->
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/MathJax.js" integrity="sha256-nvJJv9wWKEm88qvoQl9ekL2J+k/RWIsaSScxxlsrv8k=" crossorigin="anonymous"></script>
<script src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.5/config/TeX-AMS-MML_HTMLorMML.js" integrity="sha256-84DKXVJXs0/F8OTMzX4UR909+jtl4G7SPypPavF+GfA=" crossorigin="anonymous"></script>

<!--[if lt IE 9]>
<script src="https://oss.maxcdn.com/html5shiv/3.7.3/html5shiv.min.js"></script>
<script src="https://oss.maxcdn.com/respond/1.4.2/respond.min.js"></script>
<![endif]-->



  </head>

  <body data-spy="scroll" data-target="#toc">
    <div class="container template-reference-topic">
      <header>
      <div class="navbar navbar-default navbar-fixed-top" role="navigation">
  <div class="container">
    <div class="navbar-header">
      <button type="button" class="navbar-toggle collapsed" data-toggle="collapse" data-target="#navbar" aria-expanded="false">
        <span class="sr-only">Toggle navigation</span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
        <span class="icon-bar"></span>
      </button>
      <span class="navbar-brand">
        <a class="navbar-link" href="../index.html">fastTopics</a>
        <span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.5-15</span>
      </span>
    </div>

    <div id="navbar" class="navbar-collapse collapse">
      <ul class="nav navbar-nav">
        <li>
  <a href="../index.html">Home</a>
</li>
<li>
  <a href="../articles/index.html">Vignettes</a>
</li>
<li>
  <a href="../reference/index.html">Functions</a>
</li>
      </ul>
      <ul class="nav navbar-nav navbar-right">
        <li>
  <a href="https://github.com/stephenslab/fastTopics">Source</a>
</li>
      </ul>
      
    </div><!--/.nav-collapse -->
  </div><!--/.container -->
</div><!--/.navbar -->

      

      </header>

<div class="row">
  <div class="col-md-9 contents">
    <div class="page-header">
    <h1>t-SNE from Poisson NMF or Multinomial Topic Model</h1>
    <small class="dont-index">Source: <a href='https://github.com/stephenslab/fastTopics/blob/master/R/plots.R'><code>R/plots.R</code></a></small>
    <div class="hidden name"><code>tsne_from_topics.Rd</code></div>
    </div>

    <div class="ref-description">
    <p>Computes a low-dimensional nonlinear embededding of
  the data from the estimated loadings or mixture proportions using
  the t-SNE nonlinear dimensionality reduction method.</p>
    </div>

    <pre class="usage"><span class='fu'>tsne_from_topics</span>(
  <span class='no'>fit</span>,
  <span class='kw'>dims</span> <span class='kw'>=</span> <span class='fl'>2</span>,
  <span class='kw'>n</span> <span class='kw'>=</span> <span class='fl'>5000</span>,
  <span class='kw'>scaling</span> <span class='kw'>=</span> <span class='kw'>NULL</span>,
  <span class='kw'>pca</span> <span class='kw'>=</span> <span class='fl'>FALSE</span>,
  <span class='kw'>normalize</span> <span class='kw'>=</span> <span class='fl'>FALSE</span>,
  <span class='kw'>perplexity</span> <span class='kw'>=</span> <span class='fl'>100</span>,
  <span class='kw'>theta</span> <span class='kw'>=</span> <span class='fl'>0.1</span>,
  <span class='kw'>max_iter</span> <span class='kw'>=</span> <span class='fl'>1000</span>,
  <span class='kw'>eta</span> <span class='kw'>=</span> <span class='fl'>200</span>,
  <span class='kw'>check_duplicates</span> <span class='kw'>=</span> <span class='fl'>FALSE</span>,
  <span class='kw'>verbose</span> <span class='kw'>=</span> <span class='fl'>TRUE</span>,
  <span class='no'>...</span>
)</pre>

    <h2 class="hasAnchor" id="arguments"><a class="anchor" href="#arguments"></a>Arguments</h2>
    <table class="ref-arguments">
    <colgroup><col class="name" /><col class="desc" /></colgroup>
    <tr>
      <th>fit</th>
      <td><p>An object of class &#8220;poisson_nmf_fit&#8221; or
&#8220;multinom_topic_model_fit&#8221;.</p></td>
    </tr>
    <tr>
      <th>dims</th>
      <td><p>The number of dimensions in the t-SNE embedding; passed
as argument &#8220;dims&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>n</th>
      <td><p>The maximum number of rows in the loadings matrix
<code>fit$L</code> to use; when the loadings matrix has more than
<code>n</code> rows, the t-SNE embedding is computed on a random
selection of <code>n</code> rows. An upper limit on the number of rows is
used because the runtime of <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code> increases
rapidly with the number of rows in the input matrix.</p></td>
    </tr>
    <tr>
      <th>scaling</th>
      <td><p>A numeric vector of length equal to the number of
topics specifying a scaling of the columns of <code>fit$L</code>; this
re-scaling is performed prior to running t-SNE. The vector should
contain non-negative numbers only. A larger value will increase the
importance, or &#8220;weight&#8221;, of the respective topic in
computing the embedding. When <code>scaling</code> is <code>NULL</code>, no
re-scaling is performed. Note that this scaling will have no effect
if <code>normalize = TRUE</code>.</p></td>
    </tr>
    <tr>
      <th>pca</th>
      <td><p>Whether to perform a PCA processing stepe in t-SNE;
passed as argument &#8220;pca&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>normalize</th>
      <td><p>Whether to normalize the data prior to running
t-SNE; passed as argument &#8220;normalize&#8221; to
<code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>perplexity</th>
      <td><p>t-SNE perplexity parameter, passed as argument
&#8220;perplexity&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>. The perplexity
is automatically revised if it is too large; see
<code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code> for more information.</p></td>
    </tr>
    <tr>
      <th>theta</th>
      <td><p>t-SNE speed/accuracy trade-off parameter; passed as
argument &#8220;theta&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>max_iter</th>
      <td><p>Maximum number of t-SNE iterations; passed as
argument &#8220;max_iter&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>eta</th>
      <td><p>t-SNE learning rate parameter; passed as argument
&#8220;eta&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>check_duplicates</th>
      <td><p>When <code>check_duplicates = TRUE</code>, checks
whether there are duplicate rows in <code>fit$L</code>; passed as argument
&#8220;check_duplicates&#8221; to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>verbose</th>
      <td><p>If <code>verbose = TRUE</code>, progress updates are
printed; passed as argument &#8220;verbose&#8221; to
<code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    <tr>
      <th>...</th>
      <td><p>Additional arguments passed to <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>.</p></td>
    </tr>
    </table>

    <h2 class="hasAnchor" id="value"><a class="anchor" href="#value"></a>Value</h2>

    <p>A list with two list elements: <code>Y</code>, an n x d matrix
  containing the embedding <code>Y</code> returned by
  <code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code>, where n is the number of rows of the
  loadings matrix, and <code>d = dims</code>; <code>rows</code>, the rows of the
  loadings matrix included in the t-SNE embedding.</p>
    <h2 class="hasAnchor" id="details"><a class="anchor" href="#details"></a>Details</h2>

    <p>This is a lightweight interface for rapidly producing
t-SNE embeddings from matrix factorizations or multinomial topic
models; in particular, <code>tsne_from_topics</code> replaces the t-SNE
defaults with settings that are more suitable for visualizing the
structure of a matrix factorization or topic model (e.g., the PCA
step in <code>Rtsne</code> is activated by default, but disabled in
<code>tsne_from_topics</code>). See Kobak and Berens (2019) for guidance
on choosing t-SNE settings such as the "perplexity" and learning
rate (<code>eta</code>).</p>
<p>Note that since <code>tsne_plot</code> uses a <em>nonlinear</em>
transformation of the data, distances between points are less
interpretable than a linear transformation visualized using
<code>pca_plot</code> for example.</p>
    <h2 class="hasAnchor" id="references"><a class="anchor" href="#references"></a>References</h2>

    <p><NAME>. and <NAME>. (2019). The art of using t-SNE for
single-cell transcriptomics. <em>Nature Communications</em> <b>10</b>,
5416. <a href='https://doi.org/10.1038/s41467-019-13056-x'>https://doi.org/10.1038/s41467-019-13056-x</a></p>
    <h2 class="hasAnchor" id="see-also"><a class="anchor" href="#see-also"></a>See also</h2>

    <div class='dont-index'><p><code><a href='https://rdrr.io/pkg/Rtsne/man/Rtsne.html'>Rtsne</a></code></p></div>

  </div>
  <div class="col-md-3 hidden-xs hidden-sm" id="pkgdown-sidebar">
    <nav id="toc" data-toggle="toc" class="sticky-top">
      <h2 data-toc-skip>Contents</h2>
    </nav>
  </div>
</div>


      <footer>
      <div class="copyright">
  <p>Developed by <NAME>, <NAME>, <NAME>, <NAME>.</p>
</div>

<div class="pkgdown">
  <p>Site built with <a href="https://pkgdown.r-lib.org/">pkgdown</a> 1.5.1.</p>
</div>

      </footer>
   </div>

  


  </body>
</html>


