The present work shows that the  method can be applied to a sample of
neuronal recordings along two different routes: (1) apply to the sample;
or (2) apply to a larger, unsampled neuronal population from which the
sample is drawn, and then marginalize to the sample. These two routes
give inequivalent results. The second route can be further generalized to
the case where the size of the larger population is unknown. Which route
should be chosen? Some arguments are presented in favour of the second.
This work also presents and discusses probability formulae that relate
states of knowledge about a population and its samples, and that may be
useful for sampling problems in neuroscience.

§ INTRODUCTION:  AND RECORDINGS OF NEURONAL ACTIVITY

Suppose that we have recorded the firing activity of a hundred neurons,
sampled from a particular brain area. What are we to do with such data?
Gerstein, Perkel, Dayhoff <cit.> posed this question very
tersely (our emphasis):

The principal conceptual problems are (1) defining cooperativity or
functional grouping among neurons and (2) formulating
quantitative criteria for recognizing and characterizing such

These questions have a long history, of course; see for instance the 1966
review by Moore <cit.>.

The neuroscientific literature has offered several mathematical definitions
of cooperativity or functional grouping and criteria to
quantify it.

One such quantitative criterion relies on the  or relative- method
This criterion has been used in neuroscience at least since the 1990s,
applied to data recorded from brain areas as diverse as retina and motor
and it has been subjected to mathematical and conceptual scrutiny

Cooperativity can be quantified and characterized with methods in several ways. The simplest way roughly proceeds along the
following steps. Consider the recorded activity of a sample of $n$ neurons.

* The activity of each neuron, a continuous signal, is divided into $T$
time bins and binarized in intensity, and thus transformed into a
sequence of digits $0$s (inactive) and $1$s

Let the variable $\yrv_i(t)\in\set{0,1}$ denote the activity of the $i$th
sampled neuron at time bin $t$. Collectively denote the $n$ activities
with $\yr(t) \defd \bigl(\yrv_1(t),\dotsc,\yrv_n(t)\bigr)$. The
population-averaged activity at that bin is
$\yrs(t) \defd \sum_i \yrv_i(t)/n$. If we count the number of distinct
pairs of active neurons at that bin we combinatorially find
$\tbinom{n\yrs(t)}{2}\equiv n\yrs(t)\,[n\yrs(t)-1]/2$. There can be at
most $\tbinom{n}{2}$ simultaneously active pairs, so the
population-averaged pair activity is
$\av{\yr \yr}(t) \defd \tbinom{n}{2}^{-1}\tbinom{n\yrs(t)}{2}$. With some
combinatorics we see that the population-averaged activity of $m$-tuples
of neurons is

\begin{equation}
\label{eq:products_intermsof_average}
\av{\underbrace{\yr\dotsm\yr}_{\text{$m$ terms}}}(t)
= \binom{n}{m}^{-1}\binom{n\yrs(t)}{m}.
\end{equation}

For brevity let us agree to simply call activity the average
$\yrs$, pair-activity the average $\av{\yr \yr}$, and so on.

* Construct a sequence of relative- distributions for the activity
$\yrs$, using this sequence of constraints:

* the time average of the activity: $\tav{\yrs} \defd \sum_t\yrs(t)/T$;

* the time averages of the activity and of the pair-activity
$\tav{\av{\yr \yr}} \defd \sum_t\av{\yr \yr}(t)/T$;

* …
* the time averages of the activity, of the pair-activity, and so on, up
to the $k$-activity.

Call the resulting distributions $p_1(\yrs), p_2(\yrs),\dotsc,p_k(\yrs)$.
The time-bin dependence is now absent because these distributions can be
interpreted as referring to any one of the time bins $t$, or to a new
time bin (in the future or in the past) containing new data.

We also have the empirical frequency distribution of the total activity,
$f(\yrs)$, counted from the time bins.

* Now compare the distributions above with one another and with the
frequency distribution, using some probability-space distance like the
relative entropy or discrimination information
<cit.>. If we find,
say, that such distance is very high between $p_1$ and $f$, very low
between $p_2$ and $f$, and is more or less the same between all $p_m$ and
$f$ for $m \ge 2$, then we can say that there is a pairwise
cooperativity, and that any higher-order cooperativity is just a
reflection or consequence of the pairwise one. The reason is that the
information from higher-order simultaneous activities did not lead to
appreciable changes in the distribution obtained from pair activities.

The protocol above needs to be made precise by specifying various
parameters, such as the width of the time bins or the probability
distance used.

We hurry to say that the description just given is just one way to
quantify and characterize cooperativity and functional grouping, not
the only way. It can surely be criticized from many points of view.
Yet, it is quantitative and bears a more precise meaning than an undefined,
vague notion of cooperativity. Two persons who apply this
procedure to the same data will obtain the same numbers. Different
protocols can be based on the  method, for instance protocols that take
into account the activities or pair activities of specific neurons rather
than population averages, or even protocols that take into account time

The purpose of the present work is not to assess the merits of  methods
with respect to other methods. Its main purpose is to show that there is a
problem in the way the  method itself, as sketched above, is applied to
the activity of the recorded neurons. We believe that this problem is at
the root of some quirks about this method that were pointed out in the
literature <cit.>. This problems extends also to more
complex versions of the method, possibly except version that use
hidden neurons
<cit.>. The problem is
that the recorded neurons are a sample from a larger, unrecorded
population, but the  method as applied above is treating them as
isolated from the rest of the brain. Hence, the results it provides cannot
rightfully be extrapolated. We will give a mathematical proof of this. Let
us first analyse this issue in more detail.

Suppose that the neurons were recorded with electrodes covering an area of
some square millimetres <cit.>. This recording is a
sample of the activity of the neuronal population under the recording
device; which can amount to tens of thousands of neurons
<cit.>. We could even consider the recorded neurons as a sample
of a brain area more extended than the recording device.

The characterization of the cooperativity of the recorded sample would have
little meaning if we did not expect its results to generalize to a larger,
unrecorded population – at the very least that under the recording device.
In other words, we expect that the conclusions drawn with the  methods
about the sampled neurons should somehow extrapolate to unrecorded neurons
in some larger area, from which the recorded neurons were sampled. In
statistical terms we are assuming that the recorded neurons are a
representative sample of some larger neuronal population.
Probability theory tells us how to make inferences from a sample to the
larger population from which it is sampled see references

We can apply the  method to the sample, as described in the above
protocol, to generate probability distributions for the activity of the
sample. But, given that our sample is representative of a larger
population, we can also apply the  method to the larger (unrecorded)
population. The constraints are the same: the time averages of the sampled
data, since they constitute representative data about the larger population
as well. The method thus yields a probability distribution for the larger
population, and the distribution for the sample is obtained by
marginalization. The problem is that the distributions obtained from
these two applications differ. Which choice is most meaningful?

In this work we develop the second way of applying the  method, at the
level of the larger population, and show that its results differ from the
application at the sample level. We also consider the case where the size
of the larger population is unknown.

To apply the  method to the larger, unsampled population, it is
necessary to use probability relations relevant to sampling
ghoshetal1997[parts I,
VI]freedmanetal1978_r2007[ 8]gelmanetal1995_r2014[ 3]jaynes1994_r2003.
The relations we present are well-known in survey sampling and in the
pedagogic problem of drawing from an urn without replacement, yet they are
somewhat hard to find explicitly written in the neuroscientific literature.
We present and discuss them in the next section. A minor purpose of this
paper is to make these relations more widely known, because they can be
useful independently of  methods.

The notation and terminology in the present work follow ISO and ANSI standards
<cit.> but for the use of the
comma , to denote logical conjunction. Probability notation
follows Jaynes <cit.>. By probability we mean a
degree of belief which would be agreed by all rational men if
there were any rational men <cit.>.

§ PROBABILITY RELATIONS BETWEEN POPULATION AND  SAMPLE

We have already introduced the notation for the sample neurons. We
introduce an analogous notation for the $\yNv$ neurons constituting the
larger population, but using the corresponding Greek letters:
$\yRv_{\iota}(t)$ is the activity of the $\iota$th neuron at time bin $t$,
$\yRf(t) \defd \sum_{\iota} \yRv_{\iota}(t)/\yNv$ is the activity at that
bin averaged over the larger population, and so on.

The probability relations between sample and larger population are valid at
every time bin. As we mentioned above, the  distribution refers to any
time bin or to a new bin. For these reasons we will now omit the time-bin
argument (t) from our expressions.

Probabilities refer to statements about the quantities we observe. We use
the standard notation:

\begin{equation}
\label{eq:notation_statements}
\begin{aligned}
&\text{\enquote{$\yXv_{\iota} = \yRv_{\iota}$}
stands for \enquote{the activity of the $\iota$th neuron is $\yRv_{\iota}$}},
\\
&\text{\enquote{$\yXf = \yRf$}
stands for \enquote{the (population-averaged) activity of the neurons is $\yRf$}},
\\
&\text{\enquote{$\yxv_i = \yrv_i$}
stands for \enquote{the activity of the $i$th sample neuron is $\yrv_i$}},
\end{aligned}
\end{equation}

and similarly for other quantities.

If $\yHc$ denotes our state of knowledge – the evidence and assumptions
backing our probability assignments – our uncertainty about the full
activity of the larger population is expressed by the joint probability

\begin{equation}
\label{eq:joint_plaus}
\p(\yXv_1=\yRv_1, \yXv_2=\yRv_2, \dotsc, \yXv_\yNv=\yRv_\yNv \| \yHc)
\quad\text{or}\quad
\p(\yX =\yR \| \yHc), \quad \yR \in \set{0,1}^\yNv.
\end{equation}

Our uncertainty about the state of the sample is likewise expressed by

\begin{equation}
\label{eq:sample_plaus}
\p(\yxv_1=\yrv_1, \yxv_2=\yrv_2, \dotsc, \yxv_n=\yrv_n \| \yHc) \quad\text{or}\quad
\p(\yx =\yr \| \yHc), \quad \yr \in \set{0,1}^n.
\end{equation}

The theory of statistical sampling is covered in many excellent texts, for
example Ghosh  Meeden ghoshetal1997 or Freedman, Pisani, Purves [parts I, VI]freedmanetal1978_r2007; summaries can be found
in Gelman  [ 8]gelmanetal1995_r2014 and Jaynes
[ 3]jaynes1994_r2003.

We need to make an initial probability assignment for the state of the full
population before any experimental observations are made. This initial
assignment will be modified by our experimental observations, and these can
involve just a sample of the population. Our state of knowledge and initial
probability assignment should reflect that samples are somehow
representative of the whole population.

In this state of knowledge, denoted $\yH$, we know that the neurons in the
population are biologically or functionally similar, for example in
morphology or the kind of input or output they receive or give. But we are
completely ignorant about the physical details of the individual neurons.
Our ignorance is therefore symmetric under permutations of neuron
identities. This ignorance is represented by a probability distribution
that is symmetric under permutations of neuron identities; such a
distribution is usually called finitely exchangeable
ericson1969[ 1]ghoshetal1997. We stress that this
probability assignment is just an expression of the symmetry of our
ignorance about the state of the population, not an expression of
some biologic or physical symmetry or identity of the neurons.

The representation theorem for finite exchangeability states that,
in the state of knowledge $\yH$, the symmetric distribution for the full
activity is completely determined by the distribution for its

\begin{equation}
\label{eq:joint_plaus_N_homog}
\p(\yX = \yR \|  \yH) \equiv
\sum_{\yRf}\p(\yX = \yR \| \yXf = \yRf, \yH)\,
\p(\yXf=\yRf \| \yH) =
\binom{\yNv}{\yNv\yRf}^{-1} \p(\yXf=\yRf \| \yH).
\end{equation}

The equivalence on the left is just an application of the law of total
probability; the equality on the right is the statement of the theorem.
This result is intuitive: owing to symmetry, we must assign equal
probabilities to all $\tbinom{\yNv}{\yNv\yRf}$ activity vectors with
$\yNv\yRf$ active neurons; the probability of each activity vector is
therefore given by that of the average activity divided by the number of
possible vector values. Proof of this theorem and generalizations to
non-binary and continuum cases are given by de Finetti
<cit.>, Kendall <cit.>, Ericson
<cit.>, Diaconis  Freedman
<cit.>, Heath  Sudderth

Our uncertainties about the full population and the sample are connected
via the conditional probability

\begin{equation}
\label{eq:conditional_hypergeometric}
\p(\yxs = \yrs \|\yXf=\yRf, \yH)=
\binom{n}{n\yrs}\binom{\yNv-n}{\yNv \yRf-n\yrs}\binom{\yNv}{\yNv \yRf}^{-1}
\defs \ypp(\yrs\|\yRf),
\end{equation}

which is a hypergeometric distribution, typical of drawing without
replacement problems. The combinatorial proof of this expression is in
fact the same as for this class of problems
[ 3]jaynes1994_r2003[ 4.8.3]ross1976_r2010[ II.6]feller1950_r1968.

Using the conditional probability above we obtain the probability for the
activity of the sample:

\begin{gather}
\label{eq:subpop_average}
\p(\yxs=\yrs \| \yH) = \sum_{\yRf}
\p(\yxs = \yrs \|\yXf=\yRf, \yH)\,
\p(\yXf=\yRf \| \yH)
\sum_{\yRf}
\ypp(\yrs\|\yRf)\,
\p(\yXf=\yRf \| \yH).
\end{gather}

It should be proved that the probability distribution for the full activity
of the sample is also symmetric and completely determined by the
distribution of its population-averaged activity:

\begin{equation}
\label{eq:marginal}
\p(\yx = \yr \| \yH) = \binom{n}{n\yrs}^{-1} \p(\yxs=\yrs \| \yH).
\end{equation}

This is intuitively clear: our initial symmetric ignorance should also
apply to the sample. The distribution for the
sample (<ref>) indeed satisfies the same representation
theorem (<ref>) as the distribution for the full

The conditional probability
$\p(\yxs = \yrs \|\yXf=\yRf, \yH) \equiv \ypp(\yrs \|\yRf)$, besides
relating the distributions for the population and sample activities via
marginalization, also allows us to express the expectation value of any
function of the sample activity, $g(\yrs)$, in terms of the distribution
for the full population, as follows:

\begin{multline}
\label{eq:pullback_P}
\expe{g\|I}
\equiv
\sum_{\yrs} g(\yrs)\,\p(\yxs=\yrs \| \yH)
\sum_{\yrs} g(\yrs) \sum_{\yRf} \ypp(\yrs\|\yRf)\,\p(\yXf=\yRf \| \yH)
\sum_{\yRf} \biggl[ \sum_{\yrs} g(\yrs)  \ypp(\yrs\|\yRf) \biggr]\,
\p(\yXf=\yRf \| \yH),
%  =
%  \sum_{N\yRf=0}^N \sum_{n\yrs=0}^n f(\yrs)\ypp(\yrs\|\yRf)\,\p(\yXf=\yRf \| \yH)
%  \equiv  \expeb{g^*}.
\end{multline}

where the second step uses  (<ref>). The last
expression shows that the expectation of the function $g(\yrs)$ is equal to
the expectation of the function
$g^*(\yRf) \defd \sum_{\yrs} g(\yrs)\,\ypp(\yrs\|\yRf)$.

The final expression in  (<ref>) is important for our
application: the requirement that the function $g$, defined for the
sample, have a  value $c$ obtained from observed data,
translates into a linear constraint for the distribution of the full

\begin{equation}
\label{eq:constraint_extended}
c = \expe{g \| \yH} \equiv \sum_{\yRf} \biggl[ \sum_{\yrs} g(\yrs)  \ypp(\yrs\|\yRf) \biggr]\,
\p(\yXf=\yRf \| \yH).
\end{equation}

In particular, when the function $g$ is the $m$-activity of the sample,
$g(\yrs) = \sav{\yr\dotso\yr} \equiv \binom{n\yrs}{m}/\binom{n}{m}$, we

\begin{equation}
\label{eq:expe_products}
%  \label{eq:pullback_m_expectations}
%  \ypp^*\colon \sav{\underbrace{\yx \dotsm \yx}_{\text{$m$ factors}}}
%  \mapsto
% \bigl( \sav{\underbrace{\yx \dotsm \yx}_{\text{$m$ factors}}}\bigr)^* =
% \av{\underbrace{\yX \dotsm \yX}_{\text{$m$ factors}}},
% \\[2\jot]
% \expe{\sav{\underbrace{\yx \dotsm \yx}_{\text{$m$ factors}}} \| \yH} =
% \expe{\av{\underbrace{\yX \dotsm \yX}_{\text{$m$ factors}}} \| \yH}
\expe{\sav{\underbrace{\yr \dotsm \yr}_{\text{$m$ factors}}} \| \yH}
\equiv
% \\[2\jot]
% {}
% \frac{(n-m)!}{n!}
\sum_{\yrs} %m!\,
\binom{n}{m}^{-1}
\binom{n \yrs}{m}\, \p(\yxs=\yrs \| \yH)
% \sum_{N \yRf=0}^N %m!\,
% \binom{N-m}{N \yRf-m} \binom{N}{N \yRf}^{-1} \p(\yXf=\yRf \| \yH)
= \binom{\yNv}{m}^{-1}
% \frac{(N-m)!}{N!}
\sum_{\yRf} %m!\,
\binom{\yNv \yRf}{m}\, \p(\yXf=\yRf \| \yH)
\equiv    \expe{\av{\underbrace{\yR \dotsm \yR}_{\text{$m$ factors}}} \| \yH},
\end{equation}

that is, the expected values of the $m$-activities of the sample and
of the full population are equal. The proof of the middle equality uses
the expression for the $m$th factorial moment of the hypergeometric
distribution and can be found in potts1953. Similar relations
can be found for the raw moments $\expe{\yrs^m}$ and $\expe{\yRf^m}$, which
can be written in terms of the product expectations using

Thus, in a  application, when we require the expectation of the
$m$-activity of a sample to have a particular value, we are also
requiring  the expectation of the $m$-activity of the full population to
have the same value.

Log-density plot of the hypergeometric distribution
$%\p(\yxs = \yrs \|\yXf=\yRf, \yH)=
\ypp(\yrs\|\yRf) \defd  \raisebox{0pt}[0pt][1ex]{}\binom{n}{n\yrs}\binom{\yNv-n}{\yNv \yRf-n\yrs}\binom{\yNv}{\yNv \yRf}^{-1}$ for $\yNv=5000$, $n=200$. (Band artifacts may appear in the
colourbar depending on your pdf viewer.)

These expectation equalities between sample and full population should not
be surprising: we intuitively expect that the proportion of coloured
balls sampled from an urn should be roughly equal to the proportion of
coloured ball contained in the urn. The formulae in the present section
formalize and mathematically express our intuition. The hypergeometric
distribution $\ypp(\yrs\|\yRf)$ plays an important role in this
formalization. A look at its plot,  <ref>, reveals
that it is a sort of fuzzy identity transformation, or fuzzy
Kronecker delta, between the $\yRf$-space $\set{0,\dotsc,\yNv}$ and
$\yrs$-space $\set{0,\dotsc,n}$. From  (<ref>) we thus have

\begin{equation}
\label{eq:roughly_equal_nN}
\p(\yxs=a \|\yH) \approx \p(\yXf=a \|\yH),\qquad
\expeb{g(\yrs) \|\yH} \approx \expeb{g(\yRf) \|\yH},
\end{equation}

where $g$ is any smooth function defined on $\clcl{0,1}$. These approximate
equalities express the intuitive fact that our uncertainty about the
sample is representative of our uncertainty about the population and
about other samples, and vice versa. When $n=\yNv$, $\ypp(\yrs\|\yRf)$
becomes the identity matrix and the approximate equalities above become
exact – of course, since we have sampled the full population.

But the approximate equalities above may miss important features of the two
probability distributions. In the next section we will in fact emphasize
their differences. If the distribution for the population average $\yRf$ is
bimodal, for example, the bimodality can be lost in the distribution for
the sample average $\yrs$, owing to the coarsening effect of

§ MAXIMUM-ENTROPY: SAMPLE LEVEL  FULL-POPULATION LEVEL

In the previous section we have seen that observations about a sample can
be used as constraints on the distribution for the activity of the full
population. Let us use such constraints with the  method. Suppose that
we want to constrain $m$ functions of the sample activity, vectorially
written $\yf \defd (g_1,\dotsc,g_m)$, to $m$ values
$\yc \defd (c_1,\dotsc,c_m)$. These functions are typically $k$-activities
$\sav{\yr\dotso \yr}$, and the values are typically the time averages of
the observed sample, as discussed in  <ref>:
$\yc = \sum_t \yf[\yrs(t)]/T$.

Let us apply the relative- method <cit.>
directly to sampled neurons; denote this approach by $\yHb$. Then we apply
the method to the full population of neurons, most of which are unsampled;
denote this approach by $\yHa$.

Applied directly to the sampled neurons, the method yields the distribution

\begin{equation}
\label{eq:app_maxent_sample}
\p(\yxs=\yrs \|\yHb)
\binom{n}{n\yrs}
\exp[\yl\T \yf(\yrs)]
\end{equation}

where $\yk(\yl)$ is a normalization constant. The binomial in front of the
exponential appears because we must account for the multiplicity by which
the population-average activity $\yrs$ can be realized: $\yrs=0$ can be
realized in only one way (all neurons inactive), $\yrs=1/n$ can be realized
in $n$ ways (one active neuron out of $n$), and so on. This term is analogous to
the density of states in front of the Boltzmann factor in
statistical mechanics <cit.>. The $m$ Lagrange
multipliers $\yl\defd (l_1,\dotsc,l_m)$ must satisfy the $m$ constraint

\begin{equation}
\label{eq:app_maxent_sample_constraints}
c_k = \expe{g_k \|\yHb} \equiv
\frac{1}{\yk(\yl)}\sum_{\yrs}  g_k(\yrs) \binom{n}{n\yrs}
\exp[\yl\T \yf(\yrs)], \qquad k=1,\dotsc,m.
\end{equation}

Applied to the full population, using the constraint
expression (<ref>) derived in the previous section,
the method yields the distribution for the full-population activity

\begin{equation}
\label{eq:app_maxent_pop}
\p(\yXf = \yRf \| \yHa)  = \frac{1}{\yK(\yL)}\,
\binom{\yNv}{\yNv \yRf}\,\exp\Bigl[\yL\T
\sum_{\yrs} \yf(\yrs)\ypp(\yrs\|\yRf)\Bigr].
\end{equation}

The $m$ Lagrange multipliers $\yL\defd (\lambda_1,\dotsc,\lambda_m)$ must
satisfy the $m$ constraint equations

\begin{equation}
\label{eq:app_maxent_pop_constraints}
c_k = \expe{g_k \|\yHa}\equiv
\frac{1}{\yK(\yL)}\sum_{\yrs} \sum_{\yRf} g_k(\yrs) \ypp(\yrs \|\yRf)
\,
\binom{\yNv}{\yNv \yRf}\,\exp\Bigl[\yL\T
\sum_{\yrs} \yf(\yrs)\ypp(\yrs\|\yRf)\Bigr],
\qquad k=1,\dotsc,m.
\end{equation}

We obtain the distribution for the sample activity by marginalization, using

\begin{equation}
\label{eq:app_maxent_pop_marg}
\p(\yxs = \yrs \| \yHa)  = \frac{1}{\yK(\yL)}\,
\sum_{\yRf} \ypp(\yrs \|\yRf)
\,
\binom{\yNv}{\yNv \yRf}\,\exp\Bigl[\yL\T
\sum_{\yrs} \yf(\yrs)\ypp(\yrs\|\yRf)\Bigr].
\end{equation}

The distributions for the sample activity,
(<ref>) and (<ref>),
obtained with the two approaches $\yHb$ and $\yHa$, are different. From the
discussion in the previous section we expect them to be vaguely similar;
yet they cannot be exactly equal, because their equality would require the
$2m$ quantities $\yL$ and $\yl$ to satisfy the constraint equations
(<ref>) and
(<ref>), and in addition also the $n$
equations $\p(\yxs = \yrs \| \yHa) = \p(\yxs = \yrs \| \yHb)$,
$\yrs=1/n,\dotsc,1$ (one equation is taken care of by the normalization of
the distributions). We would have a set of $2m+n$ equations in $2m$

Hence, the applications of  at the sample level and at the
full-population level are inequivalent. They lead to numerically
different distributions for the sample activity $\yr$.

Linear and log-plots of $\p(\yxs = \yrs)$ constructed by  at
the population level followed by sample marginalization (blue triangles),
(<ref>), and at the sample level (red
circles),  (<ref>), with $\yNv=5000$,
$n=200$, and constraints as in  (<ref>).

Linear and log-plots of $\p(\yxs = \yrs)$ constructed by  at
the population level followed by sample marginalization (blue triangles),
(<ref>), and at the sample level (red
circles),  (<ref>), with $\yNv=5000$,
$n=200$, and constraints as in  (<ref>).

Linear and log-plots of $\p(\yxs = \yrs)$ for a sample of $n=200$
and constraints as in  (<ref>), constructed by:
red squares:  at the sample level,
blue triangles:  at the
population level,  (<ref>) with
$\yNv=10\,000$, followed by sample marginalization;
yellow circles:  at the population
level with unknown population size,
(<ref>), according to the
distribution (<ref>) for the population.

The distribution obtained at the sample level will show different
features from the one obtained at the population level, like displaced or
additional modes or particular tail behaviour. We show an example of this
discrepancy in  <ref>, for
$\yNv=10\,000$, $n=200$, and the two constraints

\begin{equation}
\label{eq:constraints}
\expe{\yrs} = 0.0478,\qquad
\expe{\yxxs} = 0.00257,
% ,\quad
% \expe{\sav{\yr\yr\yr}} = 1.48\times 10^{-4},\quad
% \expe{\sav{\yr\yr\yr\yr}} = 8.81 \times 10^{-6}.
\end{equation}

which come from the actual recording of circa 200 neurons from macaque
motor cortex <cit.>. The distribution obtained at the
population level (blue triangles) has a higher and displaced mode and a quite
different behaviour for activities around $0.5$ than the distribution
obtained at the sample level (red squares).

In our discussion we have so far assumed the size $\yNv$ of the larger
population to be known. This is rarely the case, however. We usually are
uncertain about $\yNv$ and can only guess its order of magnitude. In such a
state of knowledge $\yHd$ our ignorance about the possible value of $\yNv$
is expressed by a probability distribution $\p(\yNN=\yNv \|\yHd)=h(\yNv)$,
and the marginal distribution for the sample
activity (<ref>) is modified, by the law of total
probability, to

\begin{multline}
\label{eq:app_maxent_pop_marg_unknown_N}
\p(\yxs = \yrs \| \yHd)  =
\sum_{\yNv} \p(\yxs = \yrs \| \yNN=\yNv,\yHd) \,
\p(\yNN=\yNv \|\yHd)
\sum_\yNv \biggl\{\frac{1}{\yK(\yL_\yNv)}\,
\sum_{\yRf} \ypp_\yNv(\yrs \|\yRf)
\,
\tbinom{\yNv}{\yNv \yRf}\,\exp\Bigl[{\yL_\yNv}\T
\sum_{\yrs} \yf(\yrs)\ypp_\yNv(\yrs\|\yRf)\Bigr]\biggr\}
\,h(\yNv),
\end{multline}

where the Lagrange multipliers $\yL_{\yNv}$ and the summation range for
$\yRf$ depend on $\yNv$.

As a proof of concept,  <ref> also shows such a distribution
(yellow circles) for the same constraints as above, and a probability
distribution for $\yNv$ inspired by Jeffreys

\begin{equation}
\label{eq:pdf_popsize}
h(\yNv) \propto 1/\yNv, \qquad
\yNv \in\set{1\,000,\; 2\,000,\; \dotsc,\;10\,000}.
\end{equation}

§ DISCUSSION

The purpose of the present work was to point out and show, in a simple
set-up, that the  method can be applied to recorded neuronal data in a
way that accounts for the larger population from which the data are
sampled,  (<ref>)–(<ref>).
This application leads to results that differ from the standard application
which only considers the sample in isolation,
We gave a numerical example of this difference. We have also shown how to
extend the new application when the size of the larger population is
unknown,  (<ref>).

The latter formula, in particular, shows that the standard way of applying

implicitly assumes that no larger population exists beyond the
recorded sample of neurons. One could in fact object to the application at
the population level, and say that the traditional way of applying ,
(<ref>), yields different results because it does
not make assumptions about the size $\yNv$ of a possibly existing larger
population. Such a state of uncertainty, however, is correctly formalized
according to the laws of probability by introducing a probability
distribution for $\yNv$, and is expressed by
(<ref>). This expression cannot
generally be equal to (<ref>) unless the distribution
for $\yNv$ gives unit probability to $\yNv=n$; that is, unless the sample
is the full population, and no larger population exists.

The standard  approach therefore assumes that the recorded neurons
constitute a special subnetwork, isolated from the larger network of
neurons in which it is embedded, and which was also present under the
recording device. This assumption is unrealistic. The  approach at the
population level does not make such assumption and is therefore preferable.
It may reveal features in a data set that were unnoticed by the standard

The difference in the resulting distributions between the applications at
the sample and at the population levels appears in the use of Boltzmann
machines with hidden units <cit.>, although by a different
conceptual route. It also appears in statistical mechanics: if a system is
statistically described by a  Gibbs state, its subsystems cannot be
described by a Gibbs state <cit.>. A somewhat similar
situation also appears in the statistical description of the final state of
a non-equilibrium process starting and ending in two equilibrium states: we
can describe our knowledge about the final state either by (1) a Gibbs
distribution, calculated from the final equilibrium macrovariables, or (2)
by the distribution obtained from the Liouville evolution of the Gibbs
distribution assigned to the initial state. The two distributions differ
(even though the final physical state is obviously exactly the same
<cit.>), and the second allows us to make
sharper predictions about the final physical state thanks to our knowledge
of its preceding dynamics. In this example, though, both distributions are
usually extremely sharp and practically lead to the same predictions. In
neuroscientific applications, the difference in predictions of the sample
full-population applications can instead be very relevant.

The idea of the new application leads in fact to more questions. For

* Do the standard and new applications lead to different or contrasting
conclusions about cooperativity, when applied to real data

* How to extend the new application to the inhomogeneous case
<cit.>, in which
expectations for individual neurons or groups of neurons are constrained?

* What is the mathematical relation between the new application and
models with hidden neurons

Owing to space limitations we must leave a thorough investigation of these
questions to future work.

Finally, we would like to point out the usefulness and importance of the
probability formulae that relate our states of knowledge about a population
and its samples, presented in  <ref>. This kind of
formulae is essential in neuroscience, where we try to understand
properties of extended brain regions from partial observations. The
formulae presented here reflect a simple, symmetric state of ignorance.
More work is needed <cit.> to extend these formulae
to account for finer knowledge of the cerebral cortex and its network

§.§.§ Acknowledgements

PGLPM thanks Mari  Miri for continuous
encouragement and affection; Buster Keaton for filling life with awe and
inspiration; the developers and maintainers of , Emacs, AUC, Open
Science Framework, Python, Inkscape, Sci-Hub for making a free and
unfiltered scientific exchange possible.

§ #1

Missing 'biblatex' package
The bibliography requires the 'biblatex' package.

journaltitleJ. Neurosci.
titleCooperative firing activity in simultaneously recorded populations of neurons: detection and measurement
neuronal activity as decaying charge

journaltitleAnnu. Rev. Physiol.
titleStatistical Analysis and Functional Interpretation of Neuronal Spike Data
variable-rate Poisson process

journaltitlePhys. Rev.
note<http://bayes.wustl.edu/etj/node1.html>, see also jaynes1957b
titleInformation theory and statistical mechanics

noteRepr. in jaynes1983_r1989,  4, 39–76. <http://bayes.wustl.edu/etj/node1.html>
titleInformation theory and statistical mechanics
maximum-entropy,Bayesian probability theory,statistical fluctuations vs temporal fluctuations,phase transitions

journaltitleJ. Stat. Phys.
titleA comparison of the Shannon and Kullback information measures

Oxford University Press

noteWritten with J. Skilling. First publ. 1996
titleData Analysis: A Bayesian Tutorial

journaltitleJ. Math. Phys.
titleMaximum entropy in the problem of moments

titleMaximum entropy connections: neural networks

journaltitleBiol. Cybern.
titleDetecting higher-order interactions among the spiking events in a group of neurons

journaltitleNeural Comp.
titleThe Effects of Pair-wise and Higher-order Correlations on the Firing Rate of a Postsynaptic Neuron
maximum-entropy,computational neuroscience

journaltitleNeural Comp.
titleSynchronous Firing and Higher-Order Interactions in Neuron Pool
family=Berry II,

titleWeak pairwise correlations imply strongly correlated network states in a neural population

journaltitleJ. Neurosci.
noteSee also correction shlensetal2006_e2008
titleThe Structure of Multi-Neuron Firing Patterns in Primate Retina

Max-Planck-Institut für biologische Kybernetik

titleThe effect of pairwise neural correlations on global population statistics

journaltitlePhys. Rev. E
titleIsing model for neural data: Model quality and approximate methods for extracting functional connectivity

titleSpin glass models for a network of real neurons
maximum-entropy,Lagrange multipliers

journaltitleAdv. Neural Information Processing Systems (NIPS proceedings)
titleA joint maximum-entropy model for binary neural population patterns and continuous signals

journaltitlePhys. Rev. Lett.
titleCommon Input Explains Higher-Order Correlations and Entropy in a Simple Model of Neural Population Activity

journaltitleAdv. Neural Information Processing Systems (NIPS proceedings)
titleEmpirical models of spiking in neural populations

journaltitleProc. Natl. Acad. Sci. (USA)
titleSparse low-order interaction network underlies a highly correlated and learnable neural population code
maximum entropy

journaltitlePLoS Comput. Biol.
titleStimulus-dependent Maximum Entropy Models of Neural Population Codes

journaltitleProc. Natl. Acad. Sci. (USA)
titleThermodynamics and signatures of criticality in a network of neurons

journaltitlePhys. Rev. Lett.
titleDynamical Criticality in the Collective Activity of a Population of Retinal Neurons

journaltitleSci. Rep.
titleSimultaneous silence organizes structured higher-order interactions in neural populations

titleIsing models for networks of real neurons
maximum-entropy,Lagrange multipliers

journaltitleFront. Comput. Neurosci.
titleStatistical physics of pairwise probability models

journaltitlePLoS Comput. Biol.
titlePairwise Maximum Entropy Models for Studying Large Biological Systems: When They Can Work and When They Can't

titleWhen are microcircuits well-modeled by maximum entropy methods?
maximum-entropy,computational neuroscience

journaltitleBMC Neurosci.
numberSuppl. 1
titleWhen are microcircuits well-modeled by maximum entropy methods?
maximum-entropy,computational neuroscience

note<http://www.gatsby.ucl.ac.uk/ pel/papers/maxentbias.pdf>
titleEstimation Bias in Maximum Entropy Models

journaltitlePLoS Comput. Biol.
noteSee also the slightly different version 1605.04740
titleBistability, non-ergodicity, and inhibition in pairwise maximum-entropy models
Bayesian probability theory,computational neuroscience

journaltitleJ. Theor. Biol.
titleOutline of a theory of thought-processes and thinking machines
binarization neuronal activity

titleNeuronic Equations Revisited and Completely Solved
binarization neuronal activity

journaltitleAmerican Statistician
titleThe Kullback-Leibler distance

journaltitleJ. Stat. Phys.
titleA new theorem of information theory
probability theory,communication theory

chapter 6
titleInformation Processing in Dynamical Systems: Foundations of Harmony Theory
hidden neurons,restricted Boltzmann machines,maximum-entropy

titleCommon-input models for multiple neural spike-train data
maximum-entropy,hidden units

journaltitleJ. Phys. A
titleEffects of hidden nodes on network structure inference
hidden neurons,maximum-entropy

journaltitleJ. Phys. A
titleThe appropriateness of ignorance in the inverse kinetic Ising model
maximum-calibre,maximum-entropy,statistical mechanics,Boltzmann learning,hidden units

journaltitleJ. Neurophysiol.
titleLarge-scale, high-density (up to 512 channels) recording of local circuits in behaving animals

Cambridge University Press

titleCorticonics: Neural circuits of the cerebral cortex

seriesMonographs on statistics and applied probability
titleBayesian Methods for Finite Population Sampling
Bayesian probability theory,finite exchangeability,survey sampling

noteFirst publ. 1978
survey sampling

Boca Raton, USA

Chapman  Hall/CRC

noteFirst publ. 1995
seriesTexts in statistical science
titleBayesian Data Analysis
Bayesian probability theory,statistical models,Laplace approximation,saddle-point approximation,Monte Carlo integration

Cambridge University Press

noteEd. by G. Larry Bretthorst. First publ. 1994. <https://archive.org/details/XQUHIUXHIQUHIQXUIHX2>, <http://www-biba.inrialpes.fr/Jaynes/prob.html>, <http://omega.albany.edu:8008/JaynesBook.html>
titleProbability Theory: The Logic of Science
de Finetti's theorem,exchangeability,Bayesian probability theory

International Organization for Standardization

titleQuantities and units
standard,probability theory,terminology,mechanics,mathematics

New York

Institute of Electrical
Electronics Engineers

titleANSI/IEEE Std 260.3-1993: American National Standard: Mathematical signs and symbols for use in physical sciences and technology
Washington, D.C.

National Institute of Standards

titleGuide for the Use of the International System of Units (SI): NIST special publication 811, 1995 edition

International Organization for Standardization

titleISO 3534-1:2006: Statistics – Vocabulary and symbols – Part 1: General statistical terms and terms used in probability
standard,probability theory,terminology,repeatability,reproducibility

International Organization for Standardization

titleISO 3534-2:2006: Statistics – Vocabulary and symbols – Part 2: Applied statistics
standard,probability theory,terminology,repeatability,reproducibility

journaltitleJ. Inst. Maths. Applics
titleHow to Estimate Probabilities
Johnson-Dirichlet prior,exchangeability,frequencies of frequencies

journaltitleJ. Roy. Stat. Soc. B
note<http://www.stat.cmu.edu/ brian/905-2008/papers/Ericson-JRSSB-1969.pdf>. See also discussion sampfordetal1969
titleSubjective Bayesian Models in Sampling Finite Populations
probability theory,sufficiency,exchangeability,survey sampling

noteTransl. in definetti1972,  9, pp. 147–227
titleLa probabilità e la statistica nei rapporti con l'induzione, secondo i diversi punti di vista
Probability,Statistics and Induction: Their Relationship According to the Various Points of View,finite exchangeability

journaltitleStudia Sci. Math. Hung.
titleOn Finite and Infinite Sequences of Exchangeable Events
Bayesian probability theory,finite exchangeability,hypergeometric distribution

University of Michigan

Ann Arbor, USA

titleA Bayesian approach to two-stage sampling
partial exchangeability,de Finetti

note<http://statweb.stanford.edu/ cgates/PERSI/year.html>
titleFinite forms of de Finetti's theorem on exchangeability
Bayesian probability theory,partial exchangeability

journaltitleAnn. Prob.
titleFinite Exchangeable Sequences
de Finetti's theorem

journaltitleAmerican Statistician
titleDe Finetti's Theorem on Exchangeable Variables
Bayesian probability theory,finite exchangeability

Upper Saddle River, USA

noteFirst publ. 1976
titleA First Course in Probability
hypergeometric distribution

New York

noteFirst publ. 1950
seriesWiley series in probability and mathematical statistics
titleAn Introduction to Probability Theory and Its Applications. Vol. I

journaltitleAust. J. Phys.
titleNote on the Factorial Moments of Standard Distributions
hypergeometric distribution

New York

noteFirst publ. 1960
titleThermodynamics and an Introduction to Thermostatistics

Oxford University Press

noteFirst publ. 1939
titleTheory of Probability
Bayesian probability theory

journaltitleNeural Comp.
titleRepresentational power of restricted Boltzmann machines and deep belief networks
hidden neurons,restricted Boltzmann machines,maximum-entropy

journaltitleJ. Stat. Phys.
titleThe Restriction of the Ising Model to a Layer
non-equilibrium statistical mechanics,maximum-calibre

note<http://bayes.wustl.edu/etj/node1.html>. Extensively rewritten version of a paper first publ. 1985 in smithetal1985, pp. 377–398
titleInferential Scattering

journaltitleNat. Comm.
titleSubsampling scaling
probability theory,binomial sampling

journaltitlePhys. Rev.
note<http://bayes.wustl.edu/etj/node1.html>, see also jaynes1957
titleInformation theory and statistical mechanics. II

New York

series1962 Brandeis Summer Institute lectures in theoretical physics
titleStatistical Physics

noteEd. by R. D. Rosenkrantz. First publ. 1983
titleE. T. Jaynes: Papers on Probability, Statistics and Statistical Physics

seriesFundamental theories of physics
titleMaximum Entropy and Bayesian Methods: Laramie, Wyoming, 1990
probability theory

journaltitleJ. Neurosci.
noteSee shlensetal2006
titleCorrection, The Structure of Multi-Neuron Firing Patterns in Primate Retina

titleBrain Theory: Proceedings of the First Trieste Meeting on Brain Theory, October 1-4, 1984
neuron model
family=PDP Research Group,

Cambridge, USA

MIT Press

edition12th printing
seriesComputational models of cognition and perception
titleParallel Distributed Processing: Explorations in the Microstructure of Cognition. Vol. 1: Foundations
hidden neurons,restricted Boltzmann machines,maximum-entropy

journaltitleJ. Roy. Stat. Soc. B
note<http://www.stat.cmu.edu/ brian/905-2008/papers/Ericson-JRSSB-1969.pdf>. See ericson1969b
titleDiscussion on professor Ericson's paper
Bayesian probability theory,finite exchangeability,survey sampling

noteFirst publ. 1959
seriesCIME summer schools
titleInduzione e statistica: Lectures given at the Centro Internazionale Matematico Estivo (C.I.M.E.) held in Varenna (Como), Italy, June 1–10, 1959
finite exchangeability

titleProbability, Induction and Statistics: The art of guessing

D. Reidel

seriesFundamental theories of physics
titleMaximum-Entropy and Bayesian Methods in Inverse Problems

journaltitleJ. Roy. Stat. Soc. B
titleA Note on the Posterior Mean of a Population Mean
Bayesian probability theory,finite exchangeability,survey sampling