Bounded rationality in structured  density estimation
Authors: Tianyuan Teng
Li Kevin Wenliang
Hang Zhang
Venue: NeurIPS 2023 poster (2023)
Keywords: human representation of uncertainty; Bayesian inference; bounded rationality; inductive bias; Chinese Restaurant Process
TL;DR: We found structural inconsistencies in human density estimation and (from a large nonparametric class) we identified a bounded-rational internal model in human cognition.

Abstract: Learning to accurately represent environmental uncertainty is crucial for adaptive and optimal behaviors in various cognitive tasks. However, it remains unclear how the human brain, constrained by finite cognitive resources, constructs an internal model from an infinite space of probability distributions. In this study, we explore how these learned distributions deviate from the ground truth, resulting in observable inconsistency in a novel structured density estimation task. During each trial, human participants were asked to form and report the latent probability distribution functions underlying sequentially presented independent observations. As the number of observations increased, the reported predictive density became closer to the ground truth. Nevertheless, we observed an intriguing inconsistency in human structure estimation, specifically a large error in the number of reported clusters. Such inconsistency is invariant to the scale of the distribution and persists across stimulus modalities. We modeled uncertainty learning as approximate Bayesian inference in a nonparametric mixture prior of distributions. Human reports were best explained under resource rationality embodied in a decaying tendency towards model expansion. Our study offers insights into human cognitive processes under uncertainty and lays the groundwork for further exploration of resource-rational representations in the brain under more complex tasks.

Bounded rationality in structured density estimation
Tianyuan Teng"
Center for Life Sciences, Peking University tengtianyuan@pku. edu. cn
Li Kevin Wenliang" 3
Hang Zhang 9
Gatsby Computational Neuroscicnce Unit, University College London kevinliggatsby ucl.ucl.. ac uk
School of Psychological and Cognitive Sciences, Peking University hang. ming.zhang@pku. edu.
Abstract
Learning to accurately represent environmental uncertainty is crucial for adaptive and optimal behaviors various cognitive tasks. However, remains unclear how the human brain, constrained by finite cognitive re- sources, internalise the highly structured environmental uncortainty. In this study, explore how these learned distributions deviate from the ground truth, resulting in observable inconsistency in a novel structured density estimation task. During each trial, human participants were asked to learn and report the latent probability distribution functions underlying sequen- tially presented independent observations. As the number of observations increased, the reported predictive density became closer to the ground truth. Nevertheless, observed intriguing inconsistency in human structure estimation, specifically a large error in the number of reported clusters. Such inconsistency invariant to the scale of the distribution and persists across stimulus modalities. We modcled uncertainty lcarning as approxi- mate Bayesian inference in nonparametric mixture prior of distributions. Human reports were best explained under resource ralionality embodied in decaying lendency towards model expansion. Our study offers insights into human cognitive processes under uncertainty and lays the groundwork for further exploration of resoonce-rational representations in the brain un- der more complex tasks.
1 Introduction
We study the remarkable ability of humans to learn probabilistic distributions based on experiences in new environments, skill essential for good performance in wide range of downstream perceptual [1-4] and cognitive [5-10] tasks. This ability manifests everyday experiences, such when person learns to communicate with others of varying accents
'Equal contributions. @Co-corresponding authors. T.T. is now Max Planck Institute for Biological Cybernetics, Germany. L.K. W. is now Google DeepMind. H.Z. also affiliated with Beijing Key Laboratory of Behavior and Mental Health, PKU-IDG PKU-IDC/McGovert McGerr Institute for Brain Research, and PKU-Tsinghua Center for Life Sciences at Peking University, and Chinese Institute for Brain Research, Beijing
37th Conference on Neural Information Processing Systems (NeurIPS 2023).

---

and personalities, or discovers new social norms in foreign cultures. More impressively, humans can not only summarise key characteristics of new experiences but also describe those patterns in statistical terms, such as *Most people in X country typically talk about the weather, but some younger ones prefer to gossip instead". These observations indicate that humans can develop structured and probabilistic internal beliefs about their surroundings from experiences. To distinguish humans' internal beliefs from the modeling framework to be introduced later, refer to the former as an internal construct of the external environment. We ask: how does the human mind form internal constructs from experiences?
As experiences typically arrive sequentially in real life, building a probabilistic internal construct an online density estimation problem. An ideal solution must be both consistent in distribution estimation and invariant to the order of the experiences. Consistency requires that the learned internal construct approach the true environmental uncertainty as more data arrives, expanding in complexity when the data evidence more complicated true distribution. Invariance to experience order means no a priori preference for any particular model based on the order in which experiences arrive. For instance, if the experiences are independent samples from the same distribution, the learned distribution should not be affected by the sample order. Computationally, these two desiderata can be achieved by choosing an appropriate internal construct prior (ICP). Many instantiations of such ICPs incorporate the Chinese restaurant process (CRP) [11-13]. The CRP can recruit more resources to account for more complicated data distributions and is exchangeable so that the order of data does not affect the internal construct. In cognitive science, the CRP has been used to model category learning [14, 15, 6], classical conditioning (5, 16), among other cognitive functions [17].
However, we argue that achieving estimation consistency and order invariance is intractable given humans' limited cognitive capacity. First, buiding flexible internal construct requires an "infinite capacity" [18] that is at odds with humans' finite memory. The flexibility of the internal construct is likely bounded under resource constraints. Second, achieving order invariance online imposes significant cognitive demands: it requires memorizing all previons experiences and reallocating them into potentially different clusters for every incoming new experience. Indeed, humans do not seem to learn in an exchangeable fashion [19), which can be captured by bounded rational inference [14, 15). Resource constraints should thus impact the quality of the acquired internal construct.
To investigate how humans actually acquire and represent structured internal density model under resource constraints, we first designed a structured density estimation task wherein human participants observed sequentially presented samples drawn from a hidden Gaussian mixture with 1 clusters. Participants were asked to report the full structure of their internal construct (learned distribution) at the end of the trial, providing with high-dimensional behavioral dataset. We found, among others, that participants reported their internal constructs close to the true distribution, but they consistently reported more clusters when there were in fact fewer or even a single cluster in the true distribution.
Motivated by the resource constraint argument, we propose a flcxible ICP with variable propensity towards expansion and includes the exchangeable CRP as a special case. In addi- tion, we fit all parameters of the ICPs to the behavior dataset by novel density estimation framework (DEF), which allows us to compare between ICPs and explore other inductive biases of participants in this lask by likelihood-based mode optimization and critique. Con- sistent with our hypotheses, the proposed ICP best captures participants' behaviors when it does not preserve order-invariance and instead expands more economically when the in- ternal construct is already complex. Further, fitted ICP explains the inconsistent number of clusters in the reported constrets as having a strong prior preference for small cluster widths. In sum, our combined exporimental and modeling contributions reveal previously unknown human inductive biases in building internal constructs of new environments, and open up doors for further theoretical explorations of flexible priors in learning rich structure representations.
2

---

F
Steasw cample are
the mouse toppotte
myne repon adlonig?
Ao
NAL MM
B
C
mean (210
(om)
eant-mote 0ettrice
6.0
true dusler
<xpple
xxтрa
oventop: (9)
60
M
M
M
40
rsporte
26
M M
diuster
sample
ropo dorsily frim deraity ciseeroeboore
10
70
Figure 1: Human bchavior in a structured density estimation task. A, On each trial, human participants needed to report the latent Gaussian mixture distribution (i.c., weight, mean, and sd for each cluster) they learned from sequential samples. Panels B-F show the results of Experiment 1. B, The earth-mover distance between the reported and true distributions decreases when participants observe more samples. C, The correlation between the reported and true distributions in the first three moments. See Fig. 6E & F for detailed illustration. D, Mean number of reported clusters. Participants failed to learn the true probabilistic structure, typically reporting 2 3 clusters regardless of the true cluster number. E, Overlap ratio, which calculated as the percentage of a cluster's area being covered by any other clusters averaged across all clusters in the report. Adjacent clusters did not overlap much. F&G, Participants' reported densities versus the true distributions on example trials of Experiment 1 (F) and Experiment 2 (G). Error bars denote 1 sem over participants. Each gray line in B-D is for one participant.
Human behavioral experiments of structured density estimation
The human experiments described below were approved by the Institutional Review Board of School of Psychological and Cognitive Sciences at Peking University. All participants received monetary compensation above the local minimum wage.
2.1 Experiment 1: learning visuo-spatial distributions
We designed structured online density estimation task (Fig. 1A) to collect the internal constructs humans formed from experiences. On each trial, participant saw a sequence of red dots with horizontal positions drawn i..d. from Gaussian mixture with a varying number clusters across trials (see Appendix A.1 for the cover story and the set of distri- butions used in each experiment). At the end of the sequence, they were asked to report a Gaussian mixture model that could have generated the red dots, specifying the weight, mean, and variance of each cluster through interactive interface. As are interested in participants' innate abilities to build internal constructs, we did not provide feedback.
We first validate the quality of the reported distribution. As participants (N=21) saw more samples (70 versus 10), the earth-mover distance between the reported and true densities decreased (Fig. 1B, F(1,20) 83.21, p 0.001), which is desirable trend in density
3

---

A
stmple s2e
sampae stre
Numeric struk, small
Humeiic targe goup
sample sizo
sample siee
distunce batrewes axbucion dussers
400
diseer
Figure 2: Human learning of numeric distributions. A, Participants were randomly assigned to two groups that differed in the sd of the true distribution (scale ix vs. 5x). B, The reported cluster number. C, The distance between the adjacent clusters. Error bar denotes 1 sem over participants.
estimation. Figure 1C shows that the reported densities tracked the first three moments reasonably well (all slopes significantly greater than zero, t-tests < 0.001), although higher- order moments were harder to estimate.
A closer examination of the reported density on each trial revealed systematic inconsis- tencies in the Icarned structure. Fig. 1F shows that the reported density had more local maxima even when the true distribution only had a single Gaussian cluster. Participants tended to report around clusters, irrespective of the number of clusters in the true distri- bution (Fig. 1D). Moreover, with increasing observations, the reported number of clusters increased (1(20)=11.15 p<0..00) but did not approach the ground truth. Such inconsisten- cies in the learned number of clusters are striking and counterintuitive, especially when the true generative distribution is as simple as single Gaussian, in which case more evidence actually pushed participants further away from the truth. Another notable pattern the small overlap between the adjacent clusters (Figure 1E), which contributes to bumpiness even when the true distribution as smooth as a Gaussian. Though counterintuitive, this pattern echocs with the finding that humans near-orthogonal (i.c., non-overlapping) basis functions to represent learned visuo-motor distributions [8].
2.2 Experiments 2 and 3: distributions across domains and scales
To show that our findings in Experiment are not specific to the selected distributions, we conducted Experiment 2 (N=20) to collect human internal constructs induced by a broader set of 14 representative distributions from previous studies [20-25], with the true mumber of clusters ranging from 1 to 3. Results in Appendix A.2 show similar findings to Experiment 1: we found reasonably good density cstimation quality, and participants preferred to report 2-3 non-overlapping clusters in their internal construct, even when the true distribution is a Gaussian distribution.
Do our findings still hold when we change the domain of the distributions? In Experiment 3 (N=36), participants observed samples numerical values drawn from either Gaussian distributions or unimodal Gaussian mixtures (Fig. 2A). We divided the participants into two groups, one group observed numbers drawn from the original scale ('1X'), and the other observed samples drawn from 5 times the original scale ('5X"). After observing sequence of numbers, participants, unaware the true densitics, reported the location and weight of each latent cluster. See Appendix A.3 for details.
Consistent with Experiments and 2, we lound that participants both groups of Ex- periment 3 often reported more complex structures than the ground truth (Fig. 2B). More observations resulted in more clusters (F(1, 34) = 120.12, p < 0.001). The average reported cluster number was higher (10 samples: 2.8+0.7; 70 samples: 4.3+1.3) than Experiment 1

---

(10 samples: 2.6:0.44; 70 samples: 3.2=0.5), possibly due to differences in risuo-spatial and numeric memory capacities. Even though the distribution scales differed by a factor as large as five, both groups used a similar number of clusters to represent the numeric distribution (F(1, 34)=0.05, g=0.831). The distance between adjacent clusters differed roughly fivefold between groups, mirroring the ratio between the true scales (Fig. 2C).
Summary Through comprehensive behavioral experiments, we found that humans were able to build internal constructs of environmental uncertaintics based on online experiences, but they could not capture the number of clusters in the true distributions. Learned clusters also tended to be "orthogonal" to each other, creating unsmooth features not present in the true data distribution. It is possible that these imperfections ions related to human's limited cognitive resources, as discussed in Section 1, or specifie inductive biases in an implicit ICP, such as small cluster width. Due to the unusual complexity of the density estimation task and the nigh-dimensional human report, the detailed mechanisms are not directly obvious. We thus resort to building a mathematical model of the reported internal constructs to understand the factors contributing to these inconsistencies.
The density estimation framework
Our task explicitly queries the internal probabilistic model constructed under sequential experiences. Mapping from these experiences to distribution function is challenging and ill-posed problem. Participants need a resource-rational ICP over the space of distributions. In addition, as experimenters, we need an appropriate noise model to define a valid likelihood for maximum--lkellihon parameter fitting. To this end, we propose the density estimation framework (DEF) to model the experimental data; it consists of a rational component and an aleatoric component which we detail below. An overview is shown in Fig. 3.
3.1 The rational component
The rational component performs approximate Bayesian inference given data under an ICP, examples of which have been used extensivly in coguitive modelling [15, 16). A rational ICP for this task factorizes into prior over the latent cluster assignment of the observed dots, and a prior over the probability density function associated with each cluster. For a sequence of T dot locations Xy [x1...., .T) E RT, we denote the cluster assignments for X as 2, = [......, 21,, and each 2/ E N+ is the cluster number assigned to Fz.
The economical ICP we propose here has its prior over zy defined recursively as
Pn(2441 klzz)) = at b+or
R Eta+"
k S Kr
04 age ark, 7G := Σk1
L.kk
(1)
K2 11'
B
for E {l.... T}, where 21 = 1, n1 := Er=1 1/2 k] is the number of clusters assigned to cluster ko at time t, and K, == Σx [ 0] is the number of non-empty clusters (model size) at time The expansion rate a depends on the model size; for < 0, new clusters are less likely to be added, which implements
form of conservative expansion. The distortion rate 3 0 controls whether the different cluster sizes nG,6 are evenod-out (B < 1) or exaggerated (B > 1), form of divisive normalization. As special case, this prior recovers the conventional CRP when r = and 3 = 1, which is exchangeable and has been a popular choice for modeling flexible online learning. In contrast, the prior (1) is in general not exchangeable; see Appendix B.1.1 for examples. Also, the expected number of clusters is upper bounded by Σ, age < I-e for r > O, unlike the conventional CRP where this expectation grows as log(t). We provide justifications for the non-cschangeability of this prior at the end of this section.
The economical ICP assumees that each cluster is Gaussian distribution with mean and variance v as cluster properties. convenient prior for Gaussians is the normal-inverse--"? PR (m, v) = (m;; po, v/ An) (1;; an, oo), where the hyperparameters with subscript. 0 follow standard definitions [15]. Given XT; Bayesian inference under this TCP defined by (1) and (2) yields a posterior over Gaussian mixture parameters.
(2)
5

---

enumerate assignments
weight each by evidence
9x(2.199)
Ad 2sd) 9alKwil Az2)
create
djexp(-KK)
update
update
K
K
H
H
U
Pi+s
update
Xext
next 709
inference over an internal model prior = 6)
100
0000
000
6
0000 cluster
K
Kun
Palof
expension model modification
RIR
rational component
aleatoric component
one simulation of the density estimation framework (DEF)
simulation
Pallf (
simulation
Palpr
Σ p(()..))
X
M
BoNe perceived X
simulation M
parallel
likelihood computation
Figure 3: Schematic of the density estimation framework (DEF) and likelihood evaluation. We have omitted the subscripts & and & for clarity. The DEF comprises bounded rational component (green box) and an alcatoric component. (gold box). The former infers the cluster properties sp based on noisy perceived observations ž7 under an internal construct prior (ICP), and the latter provides a well-defined likelihood function while maintaining the dependencies between the cluster properties. Under the ICP defined in (1) and the particle- sampling algorithm in (3) and (4), the current density model either updates an existing cluster adds an additional cluster depending on the evidence of the incoming observation. The evolution of the model corresponds to a single path sampled in the tree structure shown in the left panel of the orange box. the aleatoric component, structured noise is added to ensure valid likelihood function for each simulation; see Appendix B. Each prediction yields a conditional likelihood of the reported 6p" each trial, and averaging over a large number of simulations yields the marginal likelihood. These simulations can procced in parallel (pink box), which can be accelerated by running on GPUs.
Following previous modeling work[15, 5), we assumne humans use a single particle to approx- imately perform the Bayesian updating. Specifically, the observations x; up to time and a corresponding single-particle sequence 2 (an "approximate sample" from the posterior pR (2(2(xx)) define a structured density model: a mixture of Ki Σt= 1[2] > 0] clusters, where each cluster is weighted in proportion to nf,e Σt= 1[2 ki,, and has posterior

---

over density given by the usual conjugate prior update formulae:
PR (mRs Uk Xt, zi) = -N(mii /u/ / Mnvxx ai o,
(3)
d do + mi,e.
Ni, λu +ni Ph, (Nuleo + /in // Ma
I 0000 -oooo+nnnoo ni + donie (wo - pl A0
1
di
di
where Pi, and of are, respectively, the empirical mean and the (unadjusted) empirical variance of the observations assigned to cluster & at time i, according to the particle z!. Equation (3) shows that, given zz: the sufficient statistics that determine the posterior of the k'th cluster are sp == |m, Mi d. in the sense that the conditioning variables are all functions of pi. So, we can equivalently write pR (m,, Uk Uk|X,, z)) (R Us (p))) in (3). We also regard &p (not the random variables m and in (3)) the inferred cluster properties. When the next obscrvation 24+1 arrives, it is assigned to a cluster according to the ono-step posterior
PR(2 22, Xe, Ic+1) & PR( Z))PR (2 1 B, xx) pR pR(2 (2/))// 0 241) (4)
for =(41 E [I...., Ki + 1), where pn ma(x++1/ 10 #+11 = k) is the marginal likelihood of de+1 if assigned to cluster k, given by evaluating at I141 the Student's t-density with degrees of freedom 0g, location fi, and scale (o + 1/2x,))) If 2 S K2, then the 2441 'th component is updated by T++1; if 2++1 = K2 + 1, then a new cluster created, and the internal construct expands. Figure (green box) depicts these two possibilities step at 6.
If the cluster prior (1) of the rational component implements the CRP, then by exchangeabil- ity, the true posterior over the model density is inder-invariant with respect to obscrvations. To obtain a particle from an order-invariant postcrior, one should revise the cluster assign- ments of previous x, at every time step, but this requires memorizing x, which is cognitively implausible. We argue then that keeping the cluster prior (1) exchangeable is unnecessary it forces the prior to take on cognitively implausible properly. Since (1) subsumes the CRP, we can test whether human behaviors collected in our experiments are better fit by an
prior.
3.2 The aleatoric component and maximum-likelihone model fitting
Note that each inferred particle == -ph at / = T a single simulation of the rational component a participant may possess in their mind. Our (the experimenter's) goal is to fit the rational component. to the reported cluster properties &p := ==[w,, """, aa] collecting the weights, means, and variances of the reported clusters; the reported number of clus- ters implied by sp" denoted by K", The key challenges are worth emphasizing: a) the dimensionality of the report is variable, cluster properties are unordered sets, rather than real vector; c) the cluster assignments experimenters and need to be marginalized out to obtain the
==
it depends the mumber of clusters; b) the
2y are not observed to marginal likelihood, which is intractable. We address the first two challenges by introducing aleatoric component, a structured noise model for DEF; and we deal with the last chal- lenge by an efficient implementation of the DEF that utilizes the parallel processing power of graphical processing units (GPUs).
Figure (right panel of orange box) shows illustrates the aleatoric component. Tt postulates that the participant, having inferred p', commits to the number of clusters first, during which they "slack" with small probability, reporting K clusters according to:
1, k = e, Ki E [1, (l...., Kmax) \ Ktt
PA (²)
where 0 is a slack parameter, and Kmax is the maximum number of reported clusters seen in experiment across participants. If K & K', then the participant removes splits existing clusters recursively to obtain set of slacked cluster properties ep := [n' j oi

---

(1po,K") until there are K' clusters; Appendix B.2 presents the details of how f modifies pp. Finally, the likelihood of the reported sp given slacked cluster properties are defined as P((0"(0)) == Of R + K" (infinity penalty for inferring the number of cluster wrong), and
PA = SK
Σ TESK
(w=(w";<< (1))) (o
if it = k" (6)
where Sk is the set of all permutations of K elements; (:) is Dirichlet distribution with concentration defined by n: Pp(:; is an isotropic normal with mean μ and a variance parameter; and Po (: o) is an isotropic log-normal with mean and a variance parameter. See Appendix B.2 for precise definitions of these distributions and parameters involved.
Overall, the aleatoric component the takes inferred cluster properties se from the rational component and assigns non-zero probabilities to all possible reported distributions through the slack mechanism defined by (5) and the cluster modification f. The permutation- invariant likelihood (6) regards the weight, mean, and variance vectors as independent conditional on & but retains the dependences between the clusters. To approximate the marginal likelihood, we marginalize out x' and K by averaging over a large number Monte Carlo (MC) simulations. In addition, we allow visual noise by feeding in noisy observations &, where T ~ Pa Pa( (z)) == (V(āz/ze, o.), which is also marginalized out during MC.
Altogether, the DEF estimates the likelihood of the reported sp" given xT as p(("'xx)) 2
M M E Pa " R')), Ki~ (P - (x)) & ~ (x((XT)) (7)
where M is the number of MC simulations (see Appendix B.4 for important distinction to the number of particles). We implemented this MC estimator with fully vectorized approach using PyTorch [26], which enables easy parallelization on GPUs (Fig. 3, pink box). Compared to a conventional for-loop implementation, ours yields reliable log-likelihood es- timates with orders of magnitude acceleration. This allows us to optimize the DEF by maximun-licelhhodd using any off-tho-shelve optimization method, such Nelder-Mead [27]. Further, it also allows us to critique different behavioral models (DEFs) by model comparison, using metrics such as Akaike Information Criterion (AIC) [28].
4 Experiments: fitting and comparing internal construct priors
We compare the following ICPs by fitting them in the DEFs that share the same class of alcatoric component: a) our proposed economical ICP as described in Scction 3.1; b) its special case, the exchangcable CRP-GMM; and c) a bascline batch lcarning prior that describes non-sequential cluster assignment; see Appendix B.1.3 for details. We use the Nealder-Mead oplimizer to fit the DEF parameters for each participant, restarting multiple times to avoid early convergence issues [29]. We choose large number M for MC simulations so that the estimator (7) produces small enough variance and bias tolerable for the Nelder-Mead optimizer. After fitting, we compare different. DEFs using much larger M. Details of the fitting algorithm are deferred to Appendix B.3 where we also present additional experiments validating our overall optimization scheme.
The quality of the fitted DEFs on the three experiments are shown in Fig. 4. The economical ICP produced significantly lower AICs than the exchangaalle CRP-GMM in all three ex- periments (Exp. 1: median AICs -2627.9 vs. -2583.1, < 10-4; Exp. 2: -1766.7 7 vs. -1697.9, 10-3, Exp. 3: -14.22 vs. -141.0, p < 0.01; Wilcoxon signed-rank test). We compare four aspects of the predictions in Fig. 4: the error rate in predicting the number of clusters, the negative ligeliill (NLI.) of the weights, and the expected normalized error in predict- ing mcan and oag-variance predictions using the slacked predictions given R K". The advantage of the economical ICP is not only in correctly predicting the number of clusters but also in estimating the cluster properties. As expected, the batch ICP had much worse AICs. Further, to test whether simpler models can explain the data better, we performed an ablation study based on the CRP-GMM TCP (Appendix B.5.1). The ablated ICPs ne- glect various aspects of the probabilistic structures, effectively mimicking different heuristic
8

---

error tale count
weights
eq. ertor mean
porm sqg leg-errinnce
5.0
1.4
602
Experment
0.8
1.2
100
0.6
1.0
0.4
190
1.4
103
1.2
80
0.1
1.0
-50
00
0.0
100
Exparments
60
011
thre
Figure 4: Quality of three ICPs on three experiments. Thick blue lines and whiskers are means + sem over all participants (thin lines). Red crosses are medians. Asterisk (*)
indicates significance at p<0.05 by Wilcoxon signed-rank test.
strategies in forming the internal construct. None of the ablated models outperforms the economical ICP.
Example predictions of the two sequential ICPs are shown in Fig. 5A. The DEF with the cconomical ICP captures all of the behavioral patterns exhibited by participants, including the inconsistent number of clusters, the low overlap ratio between adjacent clusters, the moments of the reported distribution (Figs. 6 to 8), and the covariance of the reported cluster propertics (Fig. 5B-D). We found largely consistent patterns in three fitted parameters across experiments (Fig. 5E 5E-G). First, in line with our resource-rational motivation, the decay rates r in (1) are greater than O, implying that model expansion is suppressed as Ot decays when the internal construct becomes more complex. Second, participants prior sd v in (2) was 20%-33 -5% of the true density's sd regardless of the true scale of the observation, with a high confidence ao (see Appendix B.5.2). This means that each cluster is kept narrow, requiring a new cluster to account for an observation far from existing clusters. It explains the bumpiness the reported internal constructs, and also suggests that participants adaptively adjusted the cluster scale in accordance with the data scale and modality. Further, the observation that participants reported around 3 clusters likely the trade-off between the opposite forces of a decaying a and a small but confident on- Third, we found the cluster distortion / (1) to be less than 1, indicating a reduced sensitivity to cluster size when forming the internal construct online.
5 Discussion
We studied density estimation by humans, sophisticated capability fundamental to va- riety of cognitive lasks but has so far been elusive in experiments and implicitly assumed in modeling. Through behavioral experiments, we discovered systematic inconsistencies in humans' internal constructs of environments (Section 2) and proposed a density estimation framework for modeling this rich behavioral dataset (Section 3). We explained the inconsis- tencies as a result of an poooonccall-expanddinn internal construct prior (ICP) over internal constructs and strong preference for narrow clusters [30). The scale-adaptive cluster prior aligns with prior rescarch on adaptive sensory and probability coding [31-33, 30, 34].
The density estimation framework (DEF) is a generic approach for modeling complex be- havioral data. Thanks to the generality, one can compare other inference algorithms based on likelihoods and can easily incorporate priors over parameters (hyperpriors). While pre- vious work kept the trainable parameters small and/or used error-based objective functions [15, 35, 4), our DEF can fit more parameters by maximum -ikelliodideohoodd marginalizing the
9

---

eccnomical
exdiongeetie
data
Gala
economicaa DEF
data
economioal DEF
clutte
DCT
4-
0.0
0.5
1.0
0.0
0.5
1,0
10
20
10
20
clsstr woight
cluster weight
cluster occentrictty
cluster accantricity
reperted clusters
catre sheeey
Redissationn
wi
Poreaurg cesavatee
neter
MAMAAA
A
MA
Eo0
Expadt
10000
Figure 5: A, Example DEF predictions for Experiment when the predicted number of clusters agree with human reports. B C, Detailed patterns of the reported cluster properties in Experiment 2, which can be successfully captured by economical DEF. In each panel, the grey heatmap denotes the distribution of responses of individual trials (collapsed across participants). The colored dots and error bars denote average response and standard error across participants in local data bins. The line denotes the regression line. B, the cluster sd increased with the cluster weight. C, The cluster sd decreased with the cluster eccentricity (distance between the cluster center and the global center). D, The cluster sd decreases with the number reported clusters. F-G, Upper panels show the effects of key parameters in the economical ICP; lower panels show the parameters fit the data of the three experiments. More visualizations of DEF predictions appear in Figs. 6 to in Appendix A
latent cluster assignments properly, thanks to our efficient implementation and hardware acceleration. With this methodology, we showed that the economical ICP provided a closer match to human behavior than the Chinese restaurant process (CRP) which has been used pervasively modeling online learning. With decaying rale of expansion, the growth of economical ICP is more in line with humans' finite memory capacity.
Our work sends important messages to the theoretical community. Dasgupta and Griffiths [32] showed that the CRP is coguitively plausible because it can be derived from preference for low entropy the cluster assignment distribution. We found that the fitted < actually increases the entropy of the prior cluster weight distribution [20, 23], suggesting that the cluster count might be the main consumption of cognitive resources. Gershman al. [16] discussed in great depth on batch versus online learning, whether the model capacily should be finite or infinite, and the hypothesis that retrospective cluster reassignment may be possible with particle representations. We provide key experimental and modeling evidence that humans may employ an online and yet finite (in expectation) model. However, with only few number of particles, the beliefs of the other possible assignments retained, so retrospective corrections may not be easily produced from these models.
not well
In sum, the cconomical DEF reproduces many human patterns, not only in the loss function optimized for (quality of capturing human structure lcarning), but also in the moments of the whole distribution and the covariance structures of cluster properties. We note the large room for improvement in predicting the cluster variance (the right column Fig. 4), which hinis al the complexity of the task and the challenge in modeling. Our work also provides comprehensive likelihood-based model comparisons on human density estimation, and paves way for future studies using high- sigh-dimensional reports in complex behavioral tasks.
10

---

Acknowledgments
This work was partly supported by the National Natural Science Foundation of China (32171095), National Science and Technology Innovation 2030 Major Program (2022ZD0204803); and funding from I.K.W. was supported by Gatsby Charitable Foundation. We thank Peter Dayan and Jian Li for helpful discussions, and five anonymous reviewers for their valuable feedback.
Center for Life Sciences H.Z.
References
[1] Marc O Ernst and Martin Banks. Humans integrate visual and haptic information in statistically optimal fashion. Nature, 2002.
[2] Yair Weiss, Eero Simoncelli, and Edward Adelson. Motion illusions as optimal percepts. Nature neuroscience, 2002.
3 Konrad P Körding and Daniel M Wolpert. Bayesian decision theory in sensorimotor control. Trends in cognitive sciences, 2006.
[4] Jake Spicer, Adam N Sanborn, and Ulrik R Beierholm. Using occam's razor and bayesian modelling to compare discrete and continuous representations in numerosity judgements. Cognitive Psychology, 2020.
5 Aaron C Courville and Nathaniel Daw. The rat as particle filter. Advances in neural information processing systems, 2007.
Thomas Griffiths, Adam N Sanborn, Kevin R Canini, and Daniel Navarro. Cat- egorization nonparametric bayesian density estimation. The probabilistic mind: Prospects for Bayesian cognitive science, 2008.
7 Nathaniel Daw and Aaron Courville. The pigeon particle filter. Advances in neural information processing systems, 2008.
8 Hang Zhang, Nathaniel D Daw, and Laurence T Maloney. Human representation of visuo-motor uncertainty as mixtures of orthogonal basis distributions. Nature Neuro- science, 2015.
99] S. J. Gershman, E. J. Horvitz, and J. B. Tenenbaum. Computational rationality: A converging paradigm for intelligence in brains, minds, and machines. Science, 2015.
[10] Kelsey Allen, Evan Shelhamer, Hanul Shin, and Joshua Tenenbaum. Infinite mixture prototypes for few-shot learning. Tn International Conference Machine Learning, 2019.
[11] David J Aldous, Mllar A Tbragimov, Jean Jacod, and David J Aldous. Eachangeability and related topics. Springer, 1985.
[12] Carl Rasmussen. The infinite gaussian mixture model. Advances in neural information processing systems, 1999.
[13] Radford M Neal. Markov chain sampling methods for dirichlet process mixture models. Journal of computational and graphical statistics, 2000.
[14] John R. Anderson. The adaptive nature of human categorization. Psychological Review, 1991.
[15] Adam N Sanborn, Thomas I Griffiths, and Daniel J Navarro. Rational approxinations to rational models: alternative algorithms for category learning. Psychological Review, 2010.
[16] Samuel Gershman, David M Blei, and Yael Niv. Context, learning, and extinction. Psychologicol Review, 2010.
11

---

[17] C. Kemp and J. B. Tenenbaum. The discovery of structural form. Proceedings of the National Academy of Sciences, 2008.
[18] Samuel J Gershman and David Blei. tutorial on bayesian nonparametric models. Journal of Mathemaaical Psychology, 2012.
[19] Bennet B Murdock Jr. The serial position effect of free recall. Journal of experimental psychology, 1962.
[20] J. Sun, J. Li, and H. Zhang. Human representation of multimodal distributions as clusters of samples. PLoS Computational Biology, 2019.
[21] Konrad Körding and Daniel M. Wolpert. The loss function of sensorimotor learning. Proceedings of the National Academy of Sciences, 2004.
[22] Konrad Körding and Daniel M Wolpert. Probabilistic Inference in Human Sensorimotor Processing. In Advances in Neural Information Processing Systems, 2003.
[23] Yeon Soon Shin and Yacl Niv. Biased evaluations emerge from inferring hidden causes. Nature Human Behaviour, 2021.
[24] T.. Acerbi, S. Vijayakumar, and D. M. Wolpert. On the origins of suboptimality in human probabilistic inference. PLoS Computational Biology, 2014.
[25] S. J. Gershman and Y. Niv. Perceptual estimation obeys Occam's razor. Frontiers in Psychology, 2013.
[26] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. Pytorch: An imperative style, high-performance deep learning library. Advances in neural infor- mation processing systems, 2019.
[27] John A Nelder and Roger Mead. A simplex method for function minimization. The compuler journal, 7(4):308 313, 1965.
[88] Hirotogu Akaike, Information theory and an extension of the maximum likelihood principle. Selected papers of hirotugu akaike, pages 199 213, 1998.
[29] William Press, Saul A Teukolsky, William T Vetterling, and Brian P Flannery. Nu- merical recipes 3rd edition: The art of scientific computing. Cambridge university press, 2007.
[30] Hang Zhang, Xiangjuan Ren, and Laurence T. Maloney. The bounded rationality of probability distortion. Proceedings of the National Academy of Sciences, 2020.
[31] Horace B. Barlow. Possible principles underlying the transformation of sensory mes- sages. Sensory communication, 1961.
[32] Ishita Dasgupta and Thomas L Griffiths. Clustering and the efficient use of cognitive resources. Journal of Mathematical Psychology. 2022.
[33] R. Polania, M. Woodford, and C. C. Ruff. Elfficinnt coding of subjective value. Nature Neuroscience, 2019.
[34] C. R. Sims. Efficient coding explains the universal law of generalization in human perception. Science, 2018
[35] Samuel J Gershman and Yael Niv. Perceptual estimation obeys occam's razor. Frontiers in psychology, 2013.
12

---

Bounded rationality in structured density estimation: Supplementary material
A Experimental details
A.1 Experiment
A.1.1 Participants
Experiment 1 recruited 21 participants (11 females, aged 18 28-25). All participants had provided informed consent before the experiment.
A.1.2 Cover story
Participants were told that they were apprentice magicians in a magical world. In this world, dangerous magic lava rocks were emitted from an unknown number of invisible volcano (cs). On cach trial, they observed past landing locations of lava rocks in specific area (on the screen), and their job was to predict the probability density of future landing locations. More specifically, they were asked to draw probability density by reporting, using click- and-drag mouse gestures, three key properties of the volcano(es), corresponding to the mean, the weight, and the standard deviation of a Gaussian component. They were told that their bonus payment depended on the accuracy of the reported predictive density.
A.1.3 Procedure and Design
On each trial, the landing positions of lava rocks were visualized red dots and sequentially presented on a black line. The number presented rocks (i.c., sample size) had two levels: 10 or 70. Each dot was presented for 166 ms, followed by a 166 ms empty screen. The landing positions were i.idd samples drawn from an unseen mixture distribution.
After observing all dots, participants needed to follow a two-stage procedure to report both (1) the predictive density of rocks' future landing positions, and (2) the underlying generative model. First, they marked the location of the volcano (es) by clicking on the black line. After marking all volcano(es). they pressed "F" to proceed to the next reporting stage. In the second stage, they chose volcano by clicking the volcano icon and then moved the mouse to report the relative density of the landing position of future lava rocks emitted from this volcano. During the report, the overall density of the rock landing position was computed and presented real-ime.
The true generative distribution set was composed of 24 distributions (Fig. 6A). The dis- tributions were created by following merge-ffom-for procedure. All distributions' overall standard deviations were about 5.3 cm. On each trial, uniform jitter ([-5.3 cm, +5.3 cml) was added the mean of the true generative distribution.
Participants completed 2 (sample size levels) x 4 (true cluster mumber levels) x 6 (distri- bution subtype) x (number of repeats) 144 trials in total. The experiment length was about 105 minutes.
A.1.4 Generation of true distribution set
The set of true distributions four-cluster Gaussian mixture distributions (bottom row in Fig. 6A) in which each com- ponent had the same SD of 0. 0.72cm of visual angle and equal weight. The three center-to- conter distances between their adjacent Gaussian components, denoted [d] dz, dg (from left to right, measured in cm), were chosen from the set of all permutations of (3.6, 4.65, 5.7)cm.
constructed following the steps below. First, created
Second, each the four-cluster Gaussian mixtures went through "merging steps", inspired by the proposal step in the dynamic clustering algorithms (e.g. Reverse-jump MCMC) in statistics ics [36). Tn each merging step, we chose two adjacent Gaussian components merge into one, with the post-merger new component having the same zeroth, first. and second
13

---

B
didritation Gauasdan
Thue
Gauanion
A
Samplescee
Sanple roll
Sur
terulesue
Thes diel Bution midure Greasiee dlasters
NM
WM
W
M
MA
M
MM
WW
MA
MM
MM
Deta
oconomical DEF
reportted mesen
Data
aconomical DEF
0.2
reported
3.0
poonomical
2.0
dAll
-10
10
10
-10
sample
sample
semple
eample
Data
economical DEF
Data
economical DEF
roported sa
cample
reported
D
Semppi sirw 10
Sample abow: 70
20
Overlep nelo (54)
15
sampla sd
sampla 5d
sample sd
sample
10
reported sherwineess
economical
reportedd skowness
economical
1ru duuster
# suue cluster
12
2
-2
sample skewness
shewness
sample skowness
Figure 6: Detailed experimental design and results in Experiment A, The distribution set. B shows the prediction of the fitted economical DEFs compared to participants' reports under different conditions. B, The relative frequeney of the reported cluster number. C, The average reported cluster number. D, The overlap ratio between the reported clusters. E, The reported moments versus the sample moments in the 10 sample-size condition, with data (left sub-panels) contrasted with model prediction (right sub-panels). Three rows are for mean, sd, and skewness. In each panel, the grey heatmap denotes the distribution of responses of individual trials (collapsed across subjects). The 5 dots and error bars denote average response and standard error across subjects in local data bins. The line denotes the regression line, with shading representing 95% confidence interval. F, The reported moments versus the sample moments in the 70 sample-size condition
moments as the combination of the two pre-merger components. The to-be-merged adjacent components were chosen in such a way that the post-merger Gaussian mixtures minimize KT- divergence KL(pprel (ippreelPpoos) (P By applying the merging step iteratively, the original four-cluster Gaussian mixtures were transformed into three-cluster, two-cluster, and finally one-cluster mixtures.
A.1.5 Results
See Fig. for the prediction of the fitted economical DEF.
A.2 Experiment
Experiment recruited participants (13 females, aged 18 25). All participants had provided informed consent before the experiment. One participant was excluded due to
14

---

Table 1: The distribution set in Experiment 2
Distribution Paper
Skewed or Symmetric
Number of clusters
Number of modes
[44]
Symmetric
Unimodal
24
Symmetric
Unimodal
25
Symmetric
Unimodal
21
Skewed
Unimodal
21
Skewed
2
Unimodal
21
Skewed
Unimodal
7
25
Skewed
2
Unimodal
23
Skewed
Bimodal
9
25
Symmetric
2
Bimodal
22
Symmetric
Bimodal
11
22
Symmetric
3
Trimodal
20
Skewed
Trimodal
13
20
Skewed
3
Trimodal
[20]
Skewed
Trimodal
This is a uniform distribution with smoothed edges.
their task performance being an outlier (measured by the earth-mover distance between the reported predictive density and the true density, exceeding 3 standard deviations, z-8core -3.3).
Experiment 1 uses 14 representative Gaussian mixture distributions selected from previ- ous studies (Table 1). The distribution set can divided into four subtypes: unimodal- symmetric, umimodal-skewed. bimodal, and trimodal. Note that the number of modes is not necessarily equal to the number of latent Gaussian clusters. For example, a unimodal skewed distribution could be mixture of two latent Gaussian clusters.
In each trial, participants observed 20 sequential samples drawn from one of the 14 distri- butions. The standard deviations of all distributions were rescaled to 4.8 cm. Each dot was presented for 166 ms, followed by a 166 ms emply screen interval. On each trial, a common uniform jitter (I-3.7 cm, +3.7 cml) was added to the means of the clusters. To balance the skewness the distribution
the experiment, we horizontally flipped the distribution in half of the trials. To explore participant's consistency in structure learning, we presented the same sample sequence for two times in seperate trials. Participants completed a total of 14 (distribution types) x 2 (flip not) (number of random sequences) 2 (number of repeats) = 112 trials. The experiment lasted about 90 minutes.
A.2.1 Results
See Fig. 7B D&G for the prediction of the economical model.
Experiment 2 contained repeated trials with identical stimuli. Using these trials, we show in Fig. 7E that in repeated trials participants reported a different number of clusters just below 50% of the time. Similarly, our model could predict with an accuracy around 0.5, very close to the participants themselves. This shows that our model can predict human report close the participants thenselves.
To illustrate the robustness of the model fitting procedure, we run model recovery exper- iments. Given randomly chosen parameters for the full model, we generate 100 sets of synthetic stimuli, reset the parameters to new random values, and then the parameters on the synthetic dataset using the procedure described in the main paper. The results show that the recovered parameters are largely consistent with the random initial values (Fig. 7F, the average correlation between the source parameters and the fitted parameters is 0.84).
A.3 Experiment
Experiment 3 recruited participants (21 females, aged 26). All participants had provided informed consent before the experiment.
15

---

A
C
20
sepporad cluater
Ovarne ralo th
Undo
2.6
modal
M
Broooi
Thenod
Unigymi
Uniaee
Tomocal
Trmmadl
Mv Mr.
M
Jnimodal-symm
Unimodal-okaa
0.8
Irocuency
0.6
report
MOSTETECall DET d312
0.-
0.2
0.0
3
reported cluster
G
Resortae aod mean
propoter ct
0.55
crick
0.00
HOSTOMECH ITE
-5
0
5
0.43
sample mean
sample mean
Namt dic aeosio
some memplee
Data
economical DEF
reported sid
10
correlation
sample sd
sample sd
0.5
reported skawnsss
economical DEF
log
60
log(nn)
logle)
sample skewness
sample skowness
Figure 7: Detailed experimental design and results for Experiment A, The distribution set. B-D&G shows the prediction of the fitted economical DEFs compared to participants' reports under different conditions. B, The relative frequency of the reported cluster number. C, The average reported cluster number. D, The overlap ratio between the reported clusters. E, The proportion that the reported mumbers of clusters are different when two trials have (1) different distributions, (2) same distribution but different samples, (3) same distribution and samples. F, The quality of model recovery, mcasured by the correlations between the source parameters and the fitted parameters. G, The reported moments versus the sample moments.
In Experiment 3, participants needed to learn the distribution of numeric values and report their belief of the distribution by entering numbers on the keyboard. On each trial, the horizontal coordinates of lava rocks were shown one-by-one on the screen, with each coordi- nate presenting for 1.5 seconds, followed by a 0.5-sccond empty screen. After observing all required to first report the number of volcanoes. Then, they were required to enter the location and the relative cruption frequency of each volcano.
coordinates, participants
The true distribution set composed of 2 unimodal distributions: one is Gaussian distribution, and the other skewed Gaussian mixture with 2 wide clusters (Fig. 2A). Participants completed 3 (distribution type) x 8 (number of repeats) 24 trials in total.
16

---

Jnimadal-Gauussaan Sample size: 10
Unimodal-skew Sample size: 10
04
frequency
report
0.2
0.0
econamical DEF
Jalmodal-GGuusiann Sample size: 70
Unimodal skew
dota
Sample size: 70
0.4
0.4
frequencc
report
0.2
0.2
0.0
0.0
reported cluster
Figure 8: The relative frequency of the reported cluster number Experiment
A.3.1 Results
See Fig. 8 for the prediction the economical model.
B
Density estimation framework details
B.1 Rational component
The rational component should take the sequential observations Xy and predict the internal construct properties s" reported by a participant. From computational perspective, any simulation-based model could be used for this component, as long as il can produce an inter- nal construct based on Xy- In this work, we adopt the rationale of Bayesian inference given a plausible generative process that could have produced the observations Xy. with approx- imations mandated by cognitive constraints. Since the behavioral report derived from subjective belief of the unobserved distribution in the environment, the generative process in the participant's mind describes subjective prior (or a construct) over distributions. As such, we call this generative process an "internal construct prior" (ICP).
The ICP in our main model inspired by the nonparametric Gaussian mixture model. A popular instantiation of the said model is defined through a CRP prior over the cluster assignments, and conjugate prior over the distribution of each cluster, abbreviated as CRP-GMM. It is appealing for our purpose because the implied generative process of xy is sequential, similar to how participants observe the xy. To make this prior more flexible, we extend the cluster assignment prior by introducing additional parameters in (1) to control the expansion decay rate (r) and count distortion rate 3. The effects of these two parameters on the prior distribution of partitions are shown in Fig. 9. While the number of clusters in the ordinary CRP grows, the decaying parameter a our extension substantially slows down the introduction of new clusters as the model size increases. As of the distortion rate 3, when B < 1, there is weaker rich-get-richer effect, giving more evenly distributed cluster sizes and higher entropy in the cluster assignment distribution, while B > produces more extreme cluster sizes and lower entropy of the cluster assignment ributinnioon.
B.1.1 Inexchangeability of the the economical ICP
The exchangeability of a cluster assignment prior is often desirable for ouline data mod- eling, because it ensures that the posterior of cluster assigument. is invariant to the order (permutation) of data. However, in this work, we are interested a prior that governs how humans construct density model online. Our discussion on cognitive const in
17

---

Section suggests that achieving order invariance requires implausible computations, so we propose relax order invariance for the purpose of modeling human cognition.
Here we show by counterexamples the proposed ICP cluster assignment prior (1) is not
exchangeable. First, for the case B and + 0. consider the partition {(21, 22): {2a)).
We have
P(21 = 12 = 1, 23 2) P = (2 = 1 = () = 21 = 1,
a0e
+ age r 2 + age
which is not the same as the probability of permuted but equivalent partition,
P( 1, 21 2, 22 2) P(=3 = 1)P(21 23 ()) 2/ 1,
ag
1 00 + oge
Then, for the case 3 * but 0, consider the partition {{21, 23, 20) (22, 21}).
We have
P(21 12 2,
1) 1 + 00
ag
2 18413
1
2-+1
2
300
2+00 00
30000
(2 + I a0
which is not the same the probability of permuted but equivalent partition,
2
200
P(21 = 23 = 1, = 1, 1,
1 00 2 00
II- a0)
Therefore, the prior defined by (1) is not exchangeable in general if B o or or #0
B.1.2 Sufficient statistics
For completeness, give the explicit expressions of the sufficient statistics in (3)
i Σ 11
PlA
Ti Σ 1 k)ar, oi.
Ti.K Σ 1| k)( - p(,x)² (8)
We assume that participants maintain these sufficient statistics when building their internal constructs, and also report them at the end of the trial, after normalizing n' to obtain the weights. Note that, ionally we denote the variance by a.
r=
decays, r=1 distortion, B= 0.5 0.5
1, 8115
0,5
number of clusters
4
2
10
20
2
10
20
2
10
20
30
2
10
20
30
2
10
20
30
time step
cluster size
20
10
2 cluster order
2
2
1 2
1 2
Figure 9: Effects of the decay rate and distortion rate B in (1) on the distribution of partitions. Heatmaps on the top row show distributions of the number of clusters as function of time steps. Brighter color means higher probability. Blue solid line indicates the mcan, blue dotted lines indicate standard deviation. The bottom row shows the distribution of cluster sizes at t = 30. Yellow bar indicates median.
18

---

B.1.3 Baseline batch ICP
As both the cconomical ICP in Section 3.1 and the proposed fitting algorithm (to appear in Appendix B.3) are new, it is important to compare the economical ICP and the exchangeable ICP, CRP-GMM, with another baseline ICP, all fit by the same algorithm. If we expect the baseline ICP to provide a worse description of human cognition than CRP-GMM on this task, then the proposed algorithm must be able to produce a less favorable model comparison result for this baseline after fitting them to human data. We choose a batch ICP as the baseline: it generates all observations as i..d. samples conditioned on cluster assignments, but, unlike the CRP, the cluster assignment prior is time-invariant and is thus exchangeable; inference over this ICP produces order-invariant posteriors. Because task has a strong sequential nature, we expect this batch ICP to be a worse descriptor of human cognition for this task.
Specifically, this batch ICP pn maintains a truncated Poisson prior distribution over Kmax submodels of mixture distributions.
(K) x Poisson((;K)) 1 5 k S K Kmaxi
(K)
(9)
10 otherwise.
The K'th submodel is a Gaussian mixture of K components. Each submodel has a sym- metric Dirichlet prior over the cluster weights, and each component has mean and variance following the conjugate Gauassian-Inverrs
distribution, as in (2). The truncated Poisson prior is heuristic choice; also tried other priors, such Categorical distribution sup- ported on K 2 and K = 3, or a truncated Poisson with mean dependent on the sample size (10 vs 70). These choices improved the ATC only insignificantly.
Inference in each submodel The following description applics to each of the Kmax submodels. Given the observed noisy data Kyy each submodel produces inferred cluster property pkx, which are the sufficient statistics of the cluster property posteriors returned by the variational-EM. algorithm (see Chapter 10.2 of [37]). We now drop the superscript i noisy
to reduce clutter, bearing in mind that all the variables are samples depending observation X and a sampled K², This algorithm alternates between the E-step: updating Dirichlet posterior over the cluster weights; and the M-step updating Gaussian- Inverse- x² posterior over the cluster means and variances. It turns out that the posterior over the cluster parameters {K depends the cluster assignments are now soft as computed by the responsibilities rt, for each F, and l E (1,..., ; K} during the variational E-step.
set of sufficient statistics similar to (8), except that
After the variational-FM procedure converges, we perform a hard cluster assignment for each observation x, by choosing the cluster with the highest responsibility,
arg max ke (E. m
Tiy.-
This ensures that the cluster weights are quantized, as is the case for the particle-based based inference methods on the other two ICPs we do not want to introduce additional effects to the variational model due to having continuous support for cluster weights Thus, at the last. time step t= T, the cluster weights are determined by пı = -'= 1[= k). The cluster mean (e and variance σg still depend on the original unquantized responsibilities and per-cluster sufficient statisties similar to (8) but with indicators 1 [|, k] replaced by responsibilities Ti, As a result, for each submodel with K clusters, the inferred cluster properties are sunmarized by ΦK := [ng, Dk, k
Submodel selection After obtaining sufficient statistics for all K models, the batch rational component selects the submodel with the largest marginal likelihood given the hard cluster assignments, computed through Student's marginal densities as in (4), minus a penalty of a weighted model size
KB = arg max
(PB(x)lςк) K
(10)
KE (e(L....mmx))
19

---

where is trainable parameter. This penalty is consistent with how AIC penalizes model complexity. Then, the submodel with KB clusters is selected, passing {ga to the aleatoric component. As such, the batch-based rational component differs from the CRP-GMM in the following ways:
1. In the batch rational component, the cluster assignments are performed by varia- tional inference rather than a particle filter;
2. the batch rational component performs model selection that penalizes large models after inferring multiple submodels, whereas CRP-GMM embeds the preference for smaller models in the cluster assignment prior.
Likelihood approximation These inferred cluster properties are passed the aleatoric component, giving (slacked) number of cluster KB and the predicted ΦB. This amounts to single simulation of the batch DEF. To compute the likelihood, the batch DEF still needs to marginalize out the visual noise and the slacked Kn in the aleatoric component. To this end, we run a large number of simulations, each with an independent draw of noisy observations Xy,, giving the predicted cluster properties Pj-- The likelihood p(p" xr) for the batch model is then approximated by (7).
The batch rational model and the aleatoric component combine to give the batch DEF, which used to benchmark the economical ICP and exchangeable ICP in the corresponding DEFs.
B.2 Aleatoric component
Here, we explain in more detail the aleatoric component described in Fig. 3. As discussed in Section 3.2, one challenge in modeling this dataset is that the dimensionality of the internal construct varies across trials. This requires that the model be able to produce variable-dimensional prodictions. In order to fit the DEFs by maximum-likclhhodd, the DEF must place nonzero probability to all possible numbers of clusters. We thus defined the distribution PA (K)K) in (5) supported on {1,..., K Kmax], and restrict the maximum number of clusters to Kaax to be the largest number of clusters ever reported by participants in an Experiment. One can also define other slack distributions over K with decaying tails to avoid an explicit upper bound.
If there is no slack, then the predicted cluster properties are as inferred. If the participant slacks with probability and commits to prediction K that does not agree with the inferred Kt from the rational component, they must modify the inferred 4p so that there are K" clusters. We assume that, during modification, the participant should keep the overall distribution roughly intact. We propose the following deterministic procodure, denoted by f fip', K'), which recursively increases or decreases the number of clusters in sp until there are Ki left.
Removing the smallest cluster. This happens whenever sp has more clusters than R We simply remove the cluster with the smallest weight. An alternative the following merging strategy: take the cluster with the smallest weight, and merge into its nearest neighbor. The merged distribution has weight equal to the sum of the weights of the clusters merged, and has mean and variances equal to the effective mean and variance of the two. More precisely, for two clusters with properties Wi, ō1] and [w2, #2, ō2] (note that we denote the variance by o), the merged cluster has properties
[w] + W2. Pm, Tm],
where pm = Wi + Wzl2 and om w1 + p + w2(62 + (2)) p.a Our results show that the removal strategy produced better AIC than the merging strategy.
Splitting the largest cluster. This happens whenever &p has fewer clusters than K', We take the cluster with the largest weight and split it into two clusters. The new clusters are centered at equal distance from and on two sides of the original cluster, and their variance is scaled version of the variance of the original cluster. Specifically, denote the properties
20

---

before fitting
after fitting
-10
14.3
est. log likelihood
-20
14.2
-30
-40
14.1
-50
14.0
-60
10a
10a
10*
10a
10a
10*
M. simulations
Figure 10: Estimated log- ligliielihood of data for one participant given different numbers of simulations. We show the mean and standard deviations of the estimate based on 30 runs for each value of M.
of this cluster to be split by [w, m, a). we the cluster propertics of the new clusters be
2 bovo, so0 2 m bova, 800
(11)
where by 0 and 0 8, free parameters to be optimised. We also tested version of this strategy where these two parameters are fixed at by gave worse AIC.
v3 and so = 1/4, but this 2
Noise processes.
The noise processes defined in (6) are given by
Pw( Dirichlet (w; c(ū 1)),
(12)
R (p(p"; II N(4e Pr: op); lem
(13)
K
P=(o";0)
II km1
log (og.Noofiss,o))
(14)
where K. implicd from G. The addition by in (12) ensures that the mode of the Dirichlet distribution equal to û. The scaling parameter e changes the confidence. The variances on and au are free parameters.
B.3 Fitting algorithm
The MC estimator for the likelihood in (7) is unbiased. However, the corresponding log- likelihood estimator, adopted in practice for numerical stability, only consistent and pro- duces nonzero bias for finite number of simulations M. This is due to Jensen's 's incquality. Suppose each simulation provides a likelioed estimate of & conditioned on some latent variables (such as xi. and Rt). Let. the true conditional likelihood be X so that ef X. then the estimated marginal lig-likelhhod is
))
E
log
log
E
-log(E[X]] Log(L[X))
(15)
M
M
Note that E[X] is the marginal likelihood. Fortunately, this bias is downwards, meaning that the expected MC estimate provides a lower bound on the true log-likelihood Suill, we must use large M during training and evaluation, as this reduces the variance of the empirical average (15) and thus lowers the bias. We show Fig. 10 the dependence of the estimated likeliood for different numbers of simulations, using data from randomly
21

---

Experiment 1
Experiment 2
Experiment 3
3.4
92
8.4
3.2
9.0
8.2
3.0
log likelhood
8.8
log Ilkeilhood
likelhood
8.0
6.6
2.6
7.8
8
2.4
7.6
22
2.0
8.0
7.4
NeelerrMead restarts
Nelder-Mead restarts
# Nelder-Mead restarts
Figure 11: Log-likelhhodd of different DEFs in ablation studies, including the economical and CRP-GMM DEFs.
picked participant in Experiment 3. When the model untrained (left panel), the bias does not completely go away even at M 10 000, but the variance is much smaller than using fewer M's. When the model is well-trained, both the bias and variance the estimated log-likelihood are substantially reduced (right pancl). Thus, the likelihood estimates become more reliable as the DEF fits the data better.
Thanks to the small estimation variance, using gradient-Tree optimization routine Nelder-Mead implemented During training, we estimate the log-likelihood using M = 10 000 parallel simulations NVIDIA GTX 1080 and A100 GPUs. The latter GPU provides likelihood estimate within 3 seconds for Experiment 2, and 5 seconds for Experiments and 3. The Nelder-Mes routine typically converges to a 0.001 relative precision on parameters within 300 iterations.
able to fit the parameters of the DEFs
the NLOpt package
We restart Nelder-Mead 10 times with parameters found from the previous optimization, as this avoids carly convergence of Nelder-Mead, To further avoid local optima, we repeat this whole procedure (with 10 restarts) 10 times with different random sceds. We can then check if our algorithm has found good solution by laking the maximum likelihood found across the 10 repeats. If this solution is reliable, then we should observe that this maximum stable across the restarts. We then take the best solution among the 10 random seeds at the last repeat as the fitted DEF parameters.
Fig. 11 shows the best (among seeds) log-likelihood across 10 restarts, averaged over all participants, for each DEF setup (see Appendix B.5.1). During the first repeat, the log- likelihood incrcases and stabilizes. Because we do not keep the best log-likllihood across the restarts, we see small fluctuations across the restarts. We see very slow increase of log- likelihood for some DEFs Experiment 3, but the relative rank of different DEFs is mostly preserved.
B.3.1 Fitting the batch DEF
For the Batch DEF, each simulation maintains Kmax submodels before the comparison in (10). Because the cluster weights posterior can be computed during variational inference, unlike in the CRP where we used particles, we do not need many simulations, and run M = 100 simulations of each of the Kmax submodels.. The output the Batch rational component is then fed into the same aleatoric component for fair comparison between the ICPs.
B.4 Particles versus MC simulations
There is an important distinction between the nunber of particles and the number of MC simulations. In this work, we assume that the participant uses single particle for infer- ring the cluster assiguments, but this is not observed from the experimenter's viewpoint. The reported internal construct may be associated with one out of many possible cluster
22

---

optimized in economical, fixed to 0 in exchangeabll
optimized in economical, fixed to 1 in exchangeable
only the baseline batch ICP only in the baseline batch ICP only the bascline batch ICP
fixed to zero, not fitted
Notes
Log-space? In default?
Y
N
Y
Y
Y
Y
N
N
Table 2: Table of all model parameters.
Y
N
Y
Y
Y
Y
Y
Y
Y
Y
Equation
(6), (12)
(6), (13)
(6), (14)
(1)
(1)
(2)
(2)
(2)
(2)
(5)
(7)
(9)
(10)
DEF component
batch rational
batch rational
batch rational
rational
rational
rational
rational
rational
rat.ioall
aleatoric
aleatoric
aleatoric
aleatoric
aleatoric
Symbol
8
A4
Ao
00
0e
Fm
K
Pseudocount of prior variance
Concentration scaling in re- porting cluster weight
Gaussian noise variance in re- porting cluster mean
Gaussian noise variance in re- porting cluster log variance
Weight model size penalty
Pseudocount of prior mean
Parameter name
Cluster count distortion
Slack on cluster count.
Visual noise variance
Cluster weight prior
Poisson prior mean
Expansion rale
Prior variance
Decay rate
Prior mean
23

---

assignments. As such, we run many simulations of the single-particle particle filter. The number of MC simulations is M, and each simulation is indexed by i. Since we stick the single-particle assumption throughout the paper, we do not need an additional index for the number of particles within each MC simulation.
Future work may consider multiple particles. In this case, each MC simulation will con- tain multiple sequences of cluster assignments. Maintaining multiple cluster assignment sequences allows re-assignment or re-weighting of the observations at each time step, achiev- ing some retrospective correction. However, is unknown how the participant makes decision on reporting the cluster properties from multiple samples, especially when the mum- ber of clusters do not agree across different particles. We thus leave multiple-particle models for future work.
B.5 Additional modeling results
After training the DEFs, we evaluate the log-likelihood using large number of simulations: 10 for the sequential DEF's in our main results, and 10 for the baseline batch DEF.
B.5.1 Ablation studies on DEF
Our initial experiment design used the exchangeable DEF (with a CRP-GMM rational com- ponent) as a reference model. We explored the model space by making small modifications to this DEF motivated by resource constraints and other heuristics. Here, we show that the two modifications in the economical DEF, namely the Decay in expansion rate a₂ and Dis- tortion in cluster count R4,4 produce reliable improvements over the reference exchangeable DEF. The DEFs resulting from other modifications cither did not produce reliable improve- ments or were insignificant compared to the reference model. All model parameters are listed in Table 2. We introduce the modifications below.
Decay: adding a model-size dependent decay rale r, in (1).
Distort: adding an exponential transformation to the size the clusters, as in (1).
Fixed Splitting: when K < K* and splitting a cluster, the parameters the splitting are fixed with b, = /3/2 and 8, = 1/4
Merge Cluster:
when K² K instcad of removing cluster, the function f(p, R) merge the smallest cluster into the cluster with closest mean. The cluster properties of the new component is computed by moment matching, as detailed in Appendix B.2.
Local MAP: instead of sampling the cluster assignment according to (4), take the cluster with largest posterior probability. This the local MAP procedure described in [15].
Constant Variance: instead of updating the variance of each cluster according to (3), the cluster mean is fixed to the trainable parameter do
Fixed Mean: instead of updating the mean of each cluster according to (3), the mean is fixed at the first obscrvation assigned to the cluster. This checks if participants are able to update the mean or simply remember the first observations assigned the clusters.
Fixed Mean Confidence: instead of fitting Ag as a parameter, we fix it at Ag = 0.01
No Visual Noise: do not add Gaussian noise to the observations.
when computing the per-step assignment posterior in (4), the like- lihood is a Student's t-distribution obtained from having a conjugate prior over the Gaussian distribution parameters. Instead of computing this likelihood using the t-distribution, this modification computes this likelihood by Gaussian with mcan Pe,k and variance dt
No Distribution Prior:
No Counting Prior: in (1), instead of using the nı, maintained for each cluster, use the average cluster size.
and Distort produced reliable improvement on ATCs compared to the CRP-GMM DEF The results of these DEFs, averaged over participants, are shown in ??. Clearly, only Decay
across all Experiments. All other modifications that deviate from rational approximation
24

---

of the Bayes rule (e.g. Fixed Mean, No Counting Prior, etc.) resulted in significantly worse fit the reported internal constructs from our participants. We also see that removing visual noise is detrimental to the quality of the fit.
The stability of the likelihood approximated from (7) is shown in Fig. 11. All DEFs fit- ting converged well as the number of restarts increases, except for a few worse DEF' in Experiment 3.
25

---

AIC
cluster
whight
3000
02
0.50
Propeoue
Dipon
Marge Clidte
Can OOW
Faod Sptting
Coustin Phior
Corstes Veesrro
No Mesat No
No entturon Proo
Propoeed
Di00t
Daca
Lood koo
Fand Spgiing
No Courtn Proe
Corf
Fass
No Moudd Nott
No nettaatzionn tror
Propoced
LOCa MSP
Mage Cluzes
Fraud Meat Corl Canics
Con Weetrro
Faod Cort
Local
Meggs
CRP
Courrio
normalzaad MSE cluster mean
normalted MSE in cluster log-variance
NLL on moan va
vaniaar varianca
015
10
0.10
Propount
Local
CRP
Food PAotSpiing
Couring Phot
Consaet Varters
to Ne'heca haree
No Dreetttt 2orr
Propout
Dipot
Decy
Looe
CaP
No Counin Pox
Redd Wean Men Conto
No Vosaaaliee
No O buton Not
Local MAP
Mouga Morgo Clustrr
CRIGMM
Spirting
Cournon Poor
Reed Mean Cald canoo
Faxd Faxt
Condlen Vettera
Phot
(a) Experiment
NIC
erter rate on dusttr count
NLL on weight
-1000
02
80
Proseted
Doo'1
Decay
Merge Cliibe
Fland Splting
Loaal M0
Food Nownn Contieers
No vouu Ne00
Propceed
Doant
Counttn Pror
CRP EAP.OM
Megge Cludee
Flonr Solltig
Riee Mexn
Corstar Varanoo
No Disiatioon Pro
No Voua Nos0
Prodosed
Draa
No Connro Proo
Solttro
MA
Find Fnod Mean Mean Confernce
Fiees Meun
Cornter Verenno
Countrg
CRP
Flend
Goratezt
No Durdtaatoon
Merg
bormalted led cluster mean
nommalzzed MSE cluster logivoranca
NLL moan vadance
0.34
.05
05
Prppord
Docay
Couring
CRP CHM
Find Speiting
Fasd Contbence
No tastt
No vssi Nose
Propasdd
Doox
Gouring Pror
CaP GVWN
Mare Chics
Food Spaiing
Local
Mear
Corstut Var
No tunin Prtr
No Woul Nose
Propoced
Decy
No Couning Frox
Iooa MMP
Fand Mean Confonnee
Toee Menn
Conater Weraro
Marp
Morg:
Faod
(b) Experiment 2
ACC
cluster count
we NLL
-100
-05
0.50
-10
Datrr
Propreed
Riss 6:01ing
Locel MAo
Conpunt Voraon Viztance
Foos Maun Cartourrs
Meepe Clasde
Do Pior
Prossced
CSP GMM
Docy
Nr Countng Pho
Local
Head Maas
No VAtu Vtoce
No Dn PHor
Prosseed
Ries Satinng
OPP CWM
Dosar
Flast Noin Confrrnnc
Seega Clidbe
Mewn
Courkng
Courtig
Local
Fiees
Merga
Condeet
nommalized MSE in cluster meen
NLL on mean verience
020
Daten
Progreed
Riess S:ktrg
GRPNEW
Lusal MMC
Corsuutt corsunt Vastanoo
Food Mour Coun Coenaann
Wer ga Chacke
Facd Moin
De Proy
Datert
Prossced
CRP-CW
Doaa
Sountin Phev
Heeld Moun
No VASLE Noor
No Sn Pfor
Looa
Foos
Mer
(c) Experiment 3
Figure 12: Model ablation comparisons for all Experiments. Lower values are better. Error bars are sems. Blue bar indicates the reference DEF with CRP-GMM as the rational component. Green bars indicate significant differences to the reference model (Wilcoxon signed-rank test, p<0..05, and grey bars indicate insignificant comparison. The Proposed (Distort + Decay) is the economical ICP.
26

---

15
economical
4
6
exchangeahle
10
4
2
0.150
0.175
0.200
0.25
0.250
0
25
75
> >00
0.0
0.2
0.4
V0
ao
Figure 13: Distribution of fitted ICP parameters using data Experiment 2.
B.5.2 Strong prior on cluster width
We found in most fitted TCPs that the prior cluster standard deviation is small (see Fig. and Fig. 13), and the confidence indicated by its pseudocount is much greater than zero. This high confidence contrasts with the low pseudocount associated with the prior cluster mean Åg. This means that participants could learn the cluster mean mostly driven by the observations while failing to adapt to the cluster uncertainty, This contributes the large number of clusters seen when there is single Gaussian in the true data distribution. We noted that Gershman and Niv [35] also used a relatively high dg = 10 parameter for their task.
References
[36] Sylvia Richardson and Peter Green On Bayesian analysis mixtures with an unknown number of components (with discussion). Journal of the Royal Statistical Society: Series (Statistical Methodology), 1997.
[37] Christopher M. Bishop. Pattern Recognition and Machine Learning, Springer, 2006.
27


Reviews:
[{"review_id": "5qoBnT0ID3", "reviewer": "Reviewer_Ds9y", "summary": "This study investigates how humans infer probability distributions from samples by combining experiments and modeling. The main contributions include a careful characterization of the behavioral tendency to overestimate the number of clusters as well as a modeling framework to identify how this behavior can arise from approximate inference. By fitting parameters from a generalized version of the Chinese restaurant process, the authors conclude that behavior can arise from a myriad of factors including strong prior expectations about the sample variance, undercounting samples, as well as a decaying tendency to form new clusters.", "questions": "Since the fitting procedure relies on approximate inference, I'd like to see a validation of the model fits to show that the method is able to recover the true parameters. Is this already published or included somewhere in the supplementary?\n\nBased on the best-fitting parameters, the authors point out a trade-off between two different effects. A decaying tendency to form new clusters (which should lead to underclustering) and a strong prior that underestimates the variance of the gaussians in the mixture (which should lead to overclustering). Could the trade-off simply be a result of over-parameterization? Again, some sort of validation of the model fitting procedure would help here.\n\nOn a related note, it wasn't clear to me how the model explains why humans report more clusters when sample size increases from 10 to 70. Shouldn't the decaying $\\alpha$ lead to fewer clusters as a function of sample size?\n\nGiven that data exchangeability is not a constraint for human inference, why not validate the model directly by running a new experiment where the sequence of samples is presented in a reverse order? I suppose the model will have concrete predictions for such an experiment.\n", "limitations": "The authors are clear about the weakness of the model in capturing aspects of the data. But they should discuss possible alternative approaches. There should be some discussion of the limitations of this paradigm and what it means for investigating structure learning in general.", "rating": 7, "confidence": 3, "soundness": 2, "presentation": 3, "contribution": 3, "strengths": "Human experiments are well-summarized and clearly document a tendency to overcluster samples. By directly asking participants to report the data-generating distribution, the approach provides a direct readout of the participant's biases that are brought to bear when solving complex probabilistic inference problems. \n\nThe paper is well-written and easy to read except for the section introducing the modeling framework where I got lost in the notations for a while. But the illustrations were very helpful.\n\nThe modeling approach is both rigorous and appropriate for answering the questions tackled by the study.", "weaknesses": "I see two main weaknesses. \n\nFirst is that the insights gained from modeling seem to be quite minimal and specific to the paradigm used here. I'd like to know whether/how the conclusions might apply more broadly to naturalistic scenarios.\n\nSecond, given that the model fits are not that great, it is quite possible that there are alternative models not considered here that explain behavior better. A discussion of alternatives is lacking.\n\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the positive feedback and constructive suggestions. \n\n***Reply to \"Implication of the results to more naturalistic scenarios\"***\n\nOur experimental design is an abstraction of many cognitive tasks that require density estimation, or finding statistical patterns from samples. A couple of motivating examples were given in lines 22-24. In addition, human’s ability to learn a probability distribution from samples is demonstrated in many prior works (Ernst & Bank 2002; Hills et al. 2002; Trommershauser et al. 2003; Kording & Wolpert 2004), but how they acquired such distributions are known.\n\nFrom a broader perspective, the current paper may help address why humans construct a biased perception of the uncertain world. For example, Ravignani et al (2016) [\\*] found that, after several rounds of laboratory music evolutions, participants’ reported distribution of rhythms became much more clustered than the original uniform distribution (Figure 2 in their paper). In other words, they found that humans exhibit an over-clustering tendency in the cultural evolution of music, which is similar to our empirical findings. Notably, Ravignani et al did not provide a cognitive model to explain why or how such over-clustering occurs in music evolution, and our DEF approach offers a promising theoretical method to analyze those results. Moreover, the overestimation of clusters we found may be related to the perception of illusory causal relationships in human superstitious thinking [\\*\\*]. Future work could consider extending our DEF to the construction of high-dimensional uncertainty representation and the emergence of causal connections among different dimensions.\n\n[\\*] Ravignani, A., Delgado, T., & Kirby, S. (2016). Musical evolution in the lab exhibits rhythmic universals. Nature Human Behaviour, 1: 0007. doi:10.1038/s41562-016-0007.\n\n[\\*\\*] Matute, H., Blanco, F., Yarritu, I., Díaz-Lago, M., Vadillo, M. A., & Barberia, I. (2015). Illusions of causality: How they bias our everyday thinking and how they could be reduced. Frontiers in Psychology, 6, 1–14.\n\n***Reply to \"Discussion of alternative models\"***\n\nPlease see the general response. Our previous statement about the model fitting performance was too pessimistic. Further results in the response PDF shows that our best model—the economical DEF captures the human data in many aspects other than cluster number and sd,  and approaches the upper limit of test-retest reliability in predicting $K$. \n\n***Reply to \"Parameter recovery\"***\n\nWe ran the suggested parameter recovery experiments and show the results in the response PDF. Given randomly chosen parameters for the full model, we generate 100 sets of synthetic stimuli, reset the parameters to new random values, and then fit the parameters on the synthetic dataset using the procedure described in the main paper. The results show that the recovered parameters are largely consistent with the random initial values (Panel K in response PDF, the average correlation between the source parameters and the fitted parameters is 0.84). We will include parameter recovery in the Appendix.\n\n***Reply to \"Overparameterization\"***\n\nThe opposite effects of decay rate and prior variance offer an explanation for the observed reported distribution on K, but these two parameters do not constitute overparameterization of the model for this given dataset. This is because the prior variance is also fit by the reported cluster variance in the dataset and is thus well-specified. There are only weak correlations between these two fitted parameters (Experiment 1: -0.12, Experiment 2: -0.05, Experiment 3: -0.05). Further, the parameter recovery experiments also indicate that these two parameters are identifiable given the data.\n\n***Reply to \"More clusters with sample size despite having a decaying alpha\"***\n\nWe clarify that the effect of a decaying alpha is to decrease the *tendency* toward creating a new cluster, rather than decreased the cluster count itself. As such, even though $\\alpha_t$ becomes small at large sample size $t$, the model size $K_t$ can still increase. \n\n***Reply to \"Additional behavioral experiments to test exchangeability\"***\n\nThe suggested experiment that tests for exchangeability by well-designed stimulus order is an interesting direction that we are currently exploring. The exact quantifiable effects of order-dependence are not directly obvious and may be intertwined with other cognitive processes, such as primacy, recency, and hindsight effects. Various hypotheses can lead to different predictions that are best tested under carefully manipulated stimulus order, opening up a large design space for exploration. We thus defer this to future work."}, {"review_id": "oX2ZwDgwq2", "reviewer": "Reviewer_Xmkb", "summary": "The authors consider the question: how do humans estimate probability distributions? They study this by performing three experiments where human subjects, where they ask participants to recover Gaussian mixture models after seeing IID samples. They find that while subjects do seem to get closer to the true mixture model as they see more samples, they consistently seem to report about three different components regardless of the ground truth number of components. They speculate that this is due to humans having limited memory.\n\nThey then propose a human model consisting of a rational component that performs approximate Bayesian inference which generalizes the standard Chinese Restaurant Process (CRP) and an _aleatoric_ component where humans pick the number of components and then merges and splits clusters until the desired number of components is reached (with some noise). They derive an efficient estimator for the likelihood under this model using Monte Carlo simulation. They then show that their model fits their experimental data much better than the standard CRP. The better performance of their model suggests that it does not seem necessary to keep the human prior exchangeable.", "questions": "1. It seems like the combination of the rational and aleatoric components are quite complicated -- have you considered simpler models? For example, instead of doing online expansion, perhaps the ICP just has a fixed prior distribution over cluster count. \n2. Similarly, given the complexity of the ICP, I'm curious how the authors interpret the fact that human subjects are able to (apparently) use it. \n3. As mentioned in the weaknesses section, I'm curious about how these results generalize to more complicated inference tasks. I understand if the experiments have not been done, but what do the authors expect the results to be if you looked at experiments with more clusters or in more realistic settings?\n4. I would've appreciated more discussion about the actual implementation of the estimator -- for example, how reliable is the log-likelihood estimate, and how sensitive is it to the number of MC simulations?\n5. Similarly, I would've appreciated more analysis of when the predictions of their proposed ICP fall short of reality (as opposed to only cases where the number of clusters agree with human report). \n6. As this is a relatively theoretically oriented paper, this is not strictly necessary. However, I wonder if the authors have seen the tendency for human subjects to pick a fixed number of clusters in other contexts, perhaps due to computational constraints. Do you expect this to have implications in practice? Are there more realistic contexts where this tendency causes systematic deviation from \"normative\" behavior?\n\nIn addition, there is the following typo in the paper:\n* Line 205-6: ... until **there** $\\hat K^i$ clusters -> ... until **there are** $\\hat K^i$ clusters", "limitations": "* I commend the authors for explicitly highlighting the relatively poor fit (albeit significantly better than prior models) of their proposed ICP on their experiments. However, I would've appreciated more discussion about why this might happen, and what next steps could be taken to address this issue.\n* I think the authors should include more discussion of the limitations and weaknesses I mentioned in the weaknesses section above. ", "rating": 8, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "* The paper is generally clearly written, and the authors do a good job of motivating their problem and describing their results. \n* The authors' structured online density estimation task is clean and straightforward. \n* The key experimental observation (that the human subjects seem to consistently estimate around 2 or 3 clusters, regardless of ground truth) is very clearly presented and clearly of great interest to this subfield of human modeling. \n* The authors two additional experiments validating the key experimental observation, with different distributions and different domains.\n* Their proposed aleatoric component follows naturally from their key observation, and cleanly captures the inconsistent number of clusters reported by human subjects. \n* The authors experimentally validate that their proposed model generally performs better than the standard CRP (as well a batch learning baseline) for all three experiments, using five different criteria: AIC, #clusters correct, negative log likelihood, the error in means, and the error in log variances.\n * The authors give the fitted parameters and discuss how to interpret them, which I found quite interesting. ", "weaknesses": "* All three experiments featured relatively simple Gaussian mixtures -- for example, experiment 1 has 1-4 clusters, experiment 2 has 1-3 clusters, and experiment 3 looks at _unimodal_ Gaussian mixtures/distributions. As a result, it seems hard to distinguish the interpretation where human subjects have limited memory and are using the proposed ICP from the interpretations where human subjects just happen to estimate 2--3 clusters across the board, or where they consistently tend to overestimate the number of clusters for other, unrelated reasons.\n* Similarly, the proposed ICP seems quite convoluted (e.g. the authors' likelihood estimator requires Monte Carlo simulation) and I'd be surprised if human subjects were actually using it.\n* Despite the complexity of the proposed ICP and the relative simplicity of the distributions being studied, the error bars for fitting human predictions are still relatively large (as shown in Figure 4), suggesting that it does not accurately capture the underlying human behavior. \n* I found figure 3 confusing; in particular, while the caption was quite helpful to my understanding, I still don't think I understand the rational or aleatoric component parts of the figure, and I did not find them helpful when I was first reading the paper. \n* Similarly, I found the discussion of how the aleatoric and rational components were composed to be initially confusing. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the positive and helpful feedback. \n\n***Reply to \"It seems hard to rule out the interpretations where human subjects just estimate 2--3 clusters across the board\"***\n\nBehaviorally, we see that the distribution of reported $K$ varies across the number of samples and distribution types (Panels B & J in PDF), so participants are not reporting $K$ randomly. Besides the batch baseline model mentioned in our general response, we also implemented another batch ICP in which the prior over $K$ is a discrete distribution over $K=2$ and $K=3$ (note that the aleatoric component still permits more values of predicted $K$). Compared to the Poisson prior in the batch baseline model, the new prior had a small but insignificant improvement (mean $\\Delta$AIC $\\approx 20.0$, $p=0.61$), still far worse than the CRP-GMM ICP. \n\n\n***Reply to \"ICP seems too complicated to be implemented by the brain\"***\n\nThis question may arise from a common confusion between the participant’s view and the experimenter’s view of the ICP. In the former, the participant always uses single particle trajectory $z_1:t$ for cluster assignment. However, these assignments are latent variables hidden away from the experimenter (us), so the experimenter needs to consider and average out all possible latent trajectories that could be generated in the participant’s mind. Only the experimenter needs to average (or marginalize) out this latent variable for modeling purposes. Note that Gershman & Niv in [25] also used a large number of simulations to marginalize out latent variables, although the likelihood function used there is more restricted.\n\n***Reply to \"Large error bars of model predictions\"***\n\nPlease see the general response, ***How well does the model produce human data? (Ds9y, Xmkb, ti72)***\n\n***Reply to \"Combination of the rational and aleatoric components\"***\n\nWe will further clarify the functions of the two components in a revision. In short, the rational component is the approximate Bayesian inverse of an ICP, implemented as a particle filter. This component produces raw predictions of the reported distribution parameters, as in most other Bayesian models for human behavior. Now, if there were only a single Gaussian in the data and in the participant’s report, then we can add simple noise models to create a valid likelihood function, such as the Dirichlet and isotropic Gaussians in Eqn (6). However, the number of clusters in the report and prediction varies in our setup, so we need to introduce a more sophisticated and structured noise model, the aleatoric component. Please see our general response for intuitions on the aleatoric component. \n\n***Reply to \"Generalisation to realistic scenarios with more clusters\"***\n\nOne key prediction of the decaying $\\alpha_t$ with $K_t$ is that humans are less likely to create a new cluster when there are already many clusters, but with overwhelming evidence from the observations, a new cluster can still be created if, for example, many clustered samples appear one next to each other and far away from existing clusters. In addition, if memorizing the clusters takes up memory capacity, then we should expect the $\\alpha_t$ to decrease faster if the dimensionality of the stimulus $x_t$ is higher since each cluster now requires multiple times more space to store. We will further verify these predictions in future work.\n\n***Reply to \"Quality of the parameter estimation procedure\"***\n\nDetailed information about the fitting algorithm is presented in Appendix B.3, where we quantify the variance and bias of the estimated log-likelihood as a function of # MC simulations in Figure 10, and we show the convergence of the parameters in Figure 11. Overall, the amount of variance and bias is very small at the optimal parameters found, although convergence is less ideal for Experiment 3. In addition, we present parameter recovery results in panel K of the response PDF. \n\n***Reply to \"Analysis of when the model fails (e.g. to predict the number of reported clusters)\"***\n\nThis is a very interesting suggestion. Note that the correlation between the statistics of reported and true distributions in Figure 1 does not exclude predictions with an unmatched number of clusters. It is also expected that the distributions would be more smooth if the predicted $K$ is less than the reported $K$, and vice versa since the cluster width is small *a priori* with high confidence. However, the behavior of the predicted cluster parameters is unspecified and unconstrained by data. Perhaps a less stringent aleatoric component may reveal more interesting predictions for an unmatched number of clusters, but any effects would be tightly linked to the inductive biases of the specific aleatoric component.\n\n***Reply to \"Tendency of humans to pick a fixed number of clusters\"*** \n\nPlease see our response to a similar question raised by **Ds9y—\"Implication of the results to more naturalistic scenarios”**.  Ravignani et al (2016) [*] found that after several rounds of laboratory music evolutions, participants’ reported distribution of rhythms became much more clustered than the original uniform distribution. According to Figure 2 in their paper, the number of clusters that emerged from participants’ reports seems to be close to 3. We also speculate that the overestimation of the number of clusters is closely related to the perception of illusory causal relationships in human superstitious thinking [**].\n\n[*] Ravignani, A., Delgado, T., & Kirby, S. (2016). Musical evolution in the lab exhibits rhythmic universals. Nature Human Behaviour, 1: 0007. doi:10.1038/s41562-016-0007.\n\n[**] Matute, H., Blanco, F., Yarritu, I., Díaz-Lago, M., Vadillo, M. A., & Barberia, I. (2015). Illusions of causality: How they bias our everyday thinking and how they could be reduced. Frontiers in Psychology, 6, 1–14.\n"}, {"review_id": "PXGj4u5vYR", "reviewer": "Reviewer_ti72", "summary": "This paper presents a visual density estimation task for human subjects based on sequential data from gaussian mixture models, complete with experimental data analyzed through a bounded-rational model of behavior. The paper reports data on the quality of the density reported by subjects, claiming that it roughly matches the true distribution's first three moments while highlighting some systematic mismatch in the number of reported mixture components. A variant of the task is performed using a different modality (numerical data instead of visual stimuli). After the task and the data is presented, the paper develops a bounded rational model for behavior. This model is composed by a generalization of the chinese restaurant process plus a sophisticated noise component. The main difference from the CRP is that this model implements in general an \"economical\" density estimation process where the likelihood of adding new clusters decreases faster than in the CRP. The behavioral model comes with an inference scheme implemented in pytorch, which allows for gpu acceleration. This model is applied to the data, showing that the \"economical\" version explains the data better than simpler alternatives (regular CRP or batch learning where cluster assignment is not sequential). This is argued to support the idea that subjects are affected by constraints on memory capacity.\n\n\n-----\n(Scores edited following the author's rebuttal)", "questions": "1. Can you provide more evidence that the quality of the density estimation reports is \"reasonable\", despite the clear systematic errors in the number of components?\n\n2. Can you clarify why the model provides experimental insight into human behavior, if it does not capture basic features of the behavior like mean and variance of the components? Moreover, how is this mismatch compatible the claim that \"The DEF with the economical ICP captures all of the behavioral patterns reported in Section 2, including the inconsistent number of clusters and the low overlap ratio between adjacent clusters\"? Figure 5A (referenced in this passage) only shows three examples; I can't find any quantitative or statistical support for this claim. Alternatively, the claims about the insights derived from the application of this model to human behavior should be significantly reduced.\n\n3. Figure captions should describe what's in the figure, not (just) provide comments. For instance, the caption of Figure 1E says \"adjacent clusters in the reported distributions did not overlap much\", but does not say what is being plotted there. How is the overlap computed? The figures would be clearer if this pattern/tendency in the captions was avoided.", "limitations": "The paper does highlight limitations in the results, although (as discussed above) other passages claim conclusions that do not seem supported by the data. Potential negative societal impact is not a concern here.", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- The paper contains a large amount of novel work, including a behavioral task, a bounded rational observer, and advanced techniques for inferring the parameters of the observer and comparing it with alternatives explanations of empirical data.\n\n- The task is ambitious in that it has a high-dimensional report.\n\n- The bounded rational observer is nontrivial but well explained, despite the limited space. From a statistical standpoint, this model and its efficient implementation are valuable.\n", "weaknesses": "1. The statements that the subjects's density estimates are reasonably good do not seem sufficiently supported by the data. Line 86: \"the reported densities tracked the ﬁrst three moments reasonably well\", but figure 1 does not show the data (the values of the true and reported moments), only the estimated correlation between them. Line 105: the density estimation quality is \"reasonably good\" in experiment 3, but figure 7 does not show any data on this - it only shows that the economical DEF reproduces certain patterns in subject behavior.\n\n2. As pointed out in the paper (lines 291-292, \"the prediction errors on the number of clusters and cluster variance are still large\"), the economical internal construct prior (ICP) model is simply not very good at describing the empirical behavioral data. This is not necessarily a problem by itself, but it makes it difficult to draw any conclusions about human behavior from the application of the model to this dataset. In particular, the claim that \"We provide key experimental and modeling evidence that humans may employ an online and yet ﬁnite (in expectation) model\" does not seem adequately supported: if the economical ICP does not describe well subject behavior, it does not matter that it does comparatively better than the CRP or batch learning.\n\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We appreciate the reviewer's suggestions on providing more details about human data and model fits, and more quantitative support for qualitative statements. These suggestions have inspired us to dig deeper into our data and model fits, which has yielded richer results (see response PDF) that we believe are more convincing. Please see detailed responses below as well as the general response to all Rs.\n\n\n***Reply to \"More evidence that the quality of the density estimation reports is reasonable\"***\n\nIn our previous manuscript, we showed that subjects’ density estimates roughly match the sample distributions in their moments (i.e., mean, sd, and skewness, Figure 1C for Experiment 1) Following the R’s suggestions, we have enhanced these plots of moments for Experiment 1 (panels D&E in response PDF) and added similar plots for Experiment 2 (panel F). In these new plots, we included three layers of data: \n1. the distribution of responses of individual trials (collapsed across subjects), plotted as the grey heatmap in the background. \n2. average responses and individual differences across subjects, plotted as the 5 black dots and their error bars (standard error).  To compute the 5 average responses, the values of the sample mean (sd, skewness) are evenly divided into 5 bins, and the trials in each bin are averaged, first separately for each subject and then across subjects.\n3. the regression lines as before.\n\nAlong with these plots of reported moments, we also present the fits of our economical DEF (sub-panels highlighted by red dots and lines). We see that the economical DEF captured not only the report of the sample mean (first row of panels D–F) but also the patterned errors in the report of sample sd and skewness (second and third rows of panels D–F). The PDF includes a few more comparisons between human data and model fits where we visualize much richer patterns in the high-dimensional human data and show how well our model fits these patterns. All of these results suggest a reasonably good fit of the economical DEF to behavioral data. \n\n***Reply to \"The estimated model does not capture data well\"***\n\nApologies for the misunderstanding caused by the original manuscript. To visualize that our model is capable of capturing various aspects of human data, we provided a more thorough illustration of model behavior in the response PDF. The fitted model behaves like the human subjects in both (1) the moments of the reported marginal density and (2) the reported cluster properties. Please see the **How well does the model produce human data? (Ds9y, Xmkb, ti72)** section in the general response for the detailed discussion.\n\n***Reply to\"What we can learn from the model if it does not capture the mean and sd.\"***\n\nThe plots in the response PDF show that the model can capture the mean, sd, and other measures in the human data. Please see the replies above as well as the general response for details.\n\n***Reply to \"Figure 5A only shows three examples and a lack of quantitative and statistical support\"***\n\nFigure 5A shows example trials, and the figures showing the overlap ratio and inconsistent report are in Figure 7 of Appendix B.4. Sorry for the slightly confusing cross-referencing. \n\n***Reply to \"Overclaims that are not well supported\"***\n\nWe will tone down the overly strong statements at the suggested places and elsewhere. We will make sure there are no other statements without data or statistical support.\n \n***Reply to \"Figure caption styles and computation of the overlap ratio\"***\n\nWe provided an intuitive visualization of the overlap ratio in the figure, but agree that we should be more explicit about the definition.  Overlap ratio is defined as the percentage of a cluster’s area being covered by any other clusters averaged across all clusters in the report.  Such computation of overlap ratio has the advantage of not being affected by the number of clusters. We also tried other measurement of the overlap ratio (e.g., the area under the second highest cluster density curve), which gave similar results. The overlap between clusters is low. We will add details of the computation of the overlap ratio to the revised text.\n\nAdding comments is common in figure captions in cognitive science papers. We will make sure that the comments are clear and faithful to the content, and that all variables/features are clearly defined in the main text. We appreciate the R’s comments.\n"}, {"review_id": "gZeWtvzFYl", "reviewer": "Reviewer_8bgA", "summary": "The authors describe a density estimation task with humans wherein participants were asked to identify parameters of an unknown distribution from presented data. While participants did a better job of estimating the overall density with more samples, the authors observed a large error in the reported number of clusters in the underlying distributions (and that increased with sample size but typically fell in the range 2-3).  The participants also tended to estimate the clusters as non-overlapping.  They modeled this process with a non-parametric mixture model.   This model consists of a “rational” component and an “aleatoric” component.  The rational component uses for cluster assignments an extension of the Dirichlet process that reduces the relative probability of adding new clusters as new points are added and modulates the relative probability of assignment to an existing cluster with a form of divisive normalization.  Unlike the regular Dirichlet process, this process is not in general exchangeable, a property the authors argue is cognitively implausible and unnecessary.  The aleatoric component adds “structured noise” in order to compute an approximation to the marginal likelihood.  They compare their model against both an exchangeable version and a baseline “batch” learning prior.", "questions": "What is the spread of data in Figure 1, particularly panels B and D?  Are participants consistent in the number of clusters assigned?  Is there a potentially simple relationship between number of clusters and sample size?  Are there simpler null models to compare against that make some relatively add hoc cluster assignments of this nature?  Maybe the baseline model behaves this way, but it is not clear to me that that is the case.  \n\nK_max from the aleatoric component is across participants, which means the posterior takes into account observations from other participants.  Is this a practical concern?  Can you replace K_max with some large but not unwieldy number and get similar results (I assume so)?\n\nThis is beyond the scope, I think, but is there a way to get running estimates of the distribution over the course of the experiment?  Your experiments don’t directly have this information, but can you get a clue to how people update models of distributions by comparing the experiments with different sample sizes?  This seems a natural question given the result that non-exchangeability seems to be important in modeling the human behavior.", "limitations": "The authors do not discuss societal impact nor directly discuss limitations.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "This is an interesting paper about the important problem of how people form internal models from observations.  The authors use an interesting approach with a non-parametric prior and provide a fitting routine that uses state-of-the-art tools for efficiency.  The experiments are clearly described and the model is (mostly) well presented.  The evidence that humans do not form models based on exchangeability is an interesting result.", "weaknesses": "Overall, the presentation is excellent, but the description and motivation of the aleatoric component is lacking, to me.  Adding some intuition about why you made the choices you made for the model would be helpful.  The “slack” probability seems particularly arbitrary to me, as does the noise added in equation 6.  Why are these necessary?  You mention in the appendix that “the DEF must place non-zero likelihood to mixture models with all possible numbers of clusters” which makes sense because many potential clusterings could result in the same reported number of clusters, but then you essentially remove from the estimation (equation 7) any sample in which the number of clusters does not match what is reported (line 208).  The “slack” probability seems to carry all of the weight of this possibility.  It is likely I’ve misunderstood something, but I think this section could be clearer.\n\nThe results from Figure 1 raise some interesting questions about the need for this particular model.  I can’t tell since they are bar plots, but it makes me wonder how consistent the cluster number estimate is and what its dependence is on the number of data points.  Moreover, can the data support a simpler psychological heuristic that assigns a number of clusters based on the number of observed samples, particularly given that subjects tend to estimate non-overlapping clusters?  It is not clear to me whether the baseline model accounts for this possibility.  I’m imagining a heuristic where people divide the screen into chunks and assign clusters that way, making finer divisions with more and more data points.  Is it clear that this is not what is happening?  Perhaps this is an ignorant question, but it seems to me that the estimation should be driven by some important visual priors.\n\n\nMinor:\n\nLine 129-130:  detailed mechanics are [not?] directly obvious…\nLine 150:  For r > 0, new clusters are… (not r< 0).\nLine 206:  There [are] \\hat{K}^i clusters…", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for helpful comments and clarifying questions. \n\nThere are many important choices in the aleatoric component, so we provide an elaborated explanation here. We will revise the manuscript and the supplemental material to add these descriptions.\n\n***Reply to \"Intuition and motivation for the aleatoric component.\"***\n\nIt is a cognitively plausible component that helps us fit the model. We regard the rational component as modeling the density estimation process while stimuli are presented, and the aleatoric component as modeling mostly the behavioral noise while participants report the cluster parameter. Please see ***Intuition and motivation for the aleatoric component*** section in our general response for a detailed reply.\n\n***Reply to \"Why is the chosen slack distribution (Eqn 5) necessary?\"***\n\nThe noise model on $K$ must assign none zero probabilities to all possible model sizes $K$ within a reasonable range, while encouraging the inferred $K^i$ from the ICP to be sampled most of the time, as it depends on the stimuli. This thus motivates the discrete K slack distribution in (5). Crucially, the simulations with incorrect $K$ are not simply removed but instead contribute zero conditional likelihood to the summation in (7). Since $M$ is fixed regardless of how many simulations yield the correct K, the ICP must produce as many as possible simulations with the correct $K$ to achieve a high likelihood. Note that the value of $\\epsilon$ does not enter the likelihood computation explicitly, but only through sampling of the slack distribution. One could also choose other conditional distributions over K for (5); however, we want to encourage the expected number of correctly predicted K to be high, rather than to encourage e.g. the expected $\\ell-1$ or $\\ell-2$ errors to be small, such as would result from a Laplace- or Gaussian-like slack distribution.\n\n***Reply to \"Why is the chosen cluster parameter distribution (6) necessary?\"***\n\nThe noise model for cluster parameters must preserve the strong dependences between \n1. the cluster parameters and reported K,\n2. between the cluster parameters themselves. \nThe latter comes from the permutation invariance of the clusters in defining the density function. These considerations motivated the seemingly complicated noise model and cluster adjustment functions in (6). Without the cluster adjustment function, the cluster parameters become decoupled from the slacked $\\hat{K}^i$. The summation over the permutations ensures that the likelihood is insensitive to the order in which the distribution parameters appear in the ordered vector representation. \n\n***Reply to \"Other heuristics where the number of clusters depends on the sample size\"***\n\nIt is possible that humans build a different prior for each sample size. However, in experimental settings where the sample size varies across trials and is thus unknown before the sequential presentation of samples (e.g. in Experiments 1 & 3), this heuristic would be incompatible with online learning. That is, agents using such heuristics need to memorize all the samples before inference, which is less cognitively plausible. For Experiment 1, we implemented a sample size-dependent prior (Poisson given each sample size) under the batch learning ICP. This gave an insignificant improvement of $\\Delta AIC=5$ (p=0.59) to batch baseline, but still far worse than the exchangeable or economical DEF. In experimental settings where the sample size is fixed (e.g., Experiment 2), this heuristic would reduce to the batch baseline model we tested in the paper and found to be inferior to the exchangeable DEF and the economical DEF. \n\n***Reply to \"Boundary-based heuristic model\"***\n\nWe thank the R for the suggested intuitive model, but we do not fully understand the suggested procedure of “making finer divisions with more and more data points”. It is possible that participants store not the clusters' means and widths, but the boundaries. However, a sensible boundary must rely on memorizing all the previous stimuli so that the divider does not cross densely populated regions. This renders this procedure less cognitively plausible. Memorizing densely populated regions is then tantamount to having a cluster-based representation, such as the GMM-based ICP we use.\n\n***Reply to \"What is the spread of data in Figure 1, particularly panels B and D?\"***\n\nTo visualize the spread, we replot Figure 1B and 1D in response PDF (Panels A & B), with individual participants’ data points denoted by grey dots linked by grey lines. \n\n***Reply to \"Are there simpler null models to compare against?\"***\n\nPlease see the general response. We also ran two additional batch ICPs using priors suggested by other reviewers. These are still significantly worse than the sequential DEFs. \n\n***Reply to \"K_max from the aleatoric component is across participants\"***\n\nThis is not a practical concern, because the amount of information provided by this $K_{max}$ is very small: only a tiny fraction of the reported $K$ reach $K_{max}$. Even if we increase $K_{max}$, then the K-slack model assigns a nonzero probability ($\\epsilon$) to a value that never appears in the dataset, so this will lower the likelihood. In this sense, $K_{max}$ is the maximum-likelihood solution of the maximum predicted $K$ for the slack model. Computationally, a $K$ too large induces a much larger number of clusters, which makes the sum in (6) over permutations impractical to compute. \n\n***Reply to \"Is there a way to get running estimates of the distribution over the course of the experiment?\"***\n\nWe considered this during our experimental design and decided not to go for this choice. Asking our participants to report on the fly may create confirmation biases or force them to commit to a distribution that is originally uncertain. Instead, we designed the 10/70 sample-size comparison without an intermediate report to test for the effect of sample size on the reported distribution.\n"}, {"review_id": "Aj80KtNywH", "reviewer": "Reviewer_ziH7", "summary": "The author built a new model to explain the human mental model of density estimation given sequentially observed data points. The model consists of a “Rational component” that generalizes the Chinese restaurant process. And an “Aleatoric component” that adds some error term. The model is fitted on real experimental data and showed superior performance compared to baseline.\n", "questions": "Experiment:\n\nI wonder how much of the over-estimation of cluster sizes is due to the small size of samples we generated? Presumably, if the ground truth number of clusters is huge, the human being will instead under-estimate the number of clusters?\n\nRational component:\n\nIt doesn’t seem too clear to me why is the process “economical”. Could you clarify this? My best understanding of the author’s argument is that when the process is exchangeable like a CRP, it involves memorizing all the points which is cognitively implausible. But I think as long as we are doing bayesian inference and getting posterior, it implies memorizing even if things aren’t exchangeable. Another interpretation is that the author is suggesting that participants are not doing a fully bayesian posterior estimation, but some step-wise procedure.\n\nI wonder if it would be useful to give some discussion of comparing economical ICP with say pitman-yor process? (e.g. some more literature review of extensions of CRP?)\n\nAleatoric component:\n\nIn general, I am a bit confused by this part. Is it just a technical point to make MLE feasible, or is it trying to model humans making mistakes?\n\nLine 194: It isn’t clear to me why “the dimensionality of the report is variable” is a challenge that needed to be solved with an aleatoric component. Could you clarify this? I thought nonparametric Bayes prior like CRP already adapts to variable dimensions well. Are you saying that the aleatoric component is definitely needed as otherwise the model cannot be fit at all? Or just that it makes things more numerically stable or realistic?\n\nFormula (5), seems a bit too simplified, is it computationally trickier to set something that gives probability based on how close K_hat is to K_i? I understand that the author suggests that we can pick whatever slack functions but I feel this particular choice somehow seems too naive.\n\nFigure 5(a): It's unclear to me the major differences between economic vs exchangeable methods?\n\nGeneral:\n\nI wonder if there should be something that models human forgetting previous examples and upweighting examples near the end. As far as I can tell that's not in the current framework.\n\nTypo:\nLine 148: n_t,k number of “clusters” assigned to cluster k. clusters should be samples?\n", "limitations": "The authors discussed about the limitations of their work in the last paragraph. I wonder if they can provide more discussions on the possible future work to address these limitations?", "rating": 5, "confidence": 2, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "The real human experimental data application looks interesting and novel. ", "weaknesses": "There are some parts of the paper that are unclear to me. See the questions section below. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for his/her insightful feedback and questions. \n\nReply to ***How much of the overestimation of cluster sizes is due to the small size of samples?***\n\nUnlikely. If the overestimation of # clusters were due to a small sample size, then we would've observed a lower # clusters in the 70-sample vs 10-sample conditions. This is not the case as we have shown in Figure 1D: the reported # clusters is higher in the 70-sample condition.\n\nReply to ***If the ground truth number of clusters is huge, the human being will instead underestimate the number of clusters?***\n\nFigure 1D suggests so: humans underestimate the number of clusters when the true number of clusters is 4. Our model also predicts so, because new clusters are less likely to form when $K_t$ is already large. \n\nReply to ***Why is the model economical?***\n\nWe thank the R for providing their interpretation of our model. The model is named economical because it is more conservative in introducing new clusters than the conventional CRP. This makes the cluster assignment prior order-dependent, as illustrated in Appendix B.1.1. We assume that this constraint is mainly due to limited memory capacity.\n\nAs the R pointed out, we do not perform full Bayesian posterior estimation. However, this feature is shared between the economical ICP and exchangeable ICP, so this is not the main reason for calling the model we proposed as being economical. \n\nReply to ***Comparing economical ICP with other extensions of CRP***\n\nThe suggested Pitman-Yor process (PYP) has a discounting parameter to control how fast new clusters are added. At each time $t$,\n$$\nP(z_{t+1}=k|\\mathbb{z}_{t}) := \\frac{\\alpha + n_k b}{\\alpha+t}, \\quad  k = K_t + 1\n$$\n\nor \n\n$$\nP(z_{t+1}=k|\\mathbb{z}_{t}) := \\frac{t - b}{\\alpha+t}, \\quad  k \\le K_t\n$$\n \nTo slow down cluster creation, we must make $b<0$. As $K_t$ grows, the probability of expansion becomes negative, posing a more challenging constrained optimization problem. The economical ICP modulates only the expansion rate through a monotonic decreasing function of $K_t$, so it is more flexible, interpretable, and easier to optimize.\n\nWe are also aware of other priors, such as the sticky CRP [*] used for memory updating [**]. However, we did not find this parametrization to improve the model fit to data in our earlier experiments. We will discuss the PYP, the sticky CRP, and other CRP-like variants in a revision. \n\n[*] Fox E B, Sudderth E B, Jordan M I, et al. A sticky HDP-HMM with application to speaker diarization. The Annals of Applied Statistics, 2011\n[**] Gershman S J, Radulescu A, Norman K A, et al. Statistical computations underlying the dynamics of memory updating. PLoS computational biology, 2014\n\nReply to ***Purpose of the aleatoric component***\n\nThe aleatoric component can be regarded as modeling human mistakes. Note that cognitive and behavioral noises are observed in many cognitive science studies, so the including the noise component is a standard procedure in modeling. See the ***Intuition and motivation for the aleatoric component (8bgA, ziH7)*** section in our general response for a detailed reply. We are happy to provide further clarification if the R could elaborate on their confusion.\n\nReply to ***Variable dimensionality***\n\nWe clarify that the challenge is not how to produce variable-dimensional predictions, but rather how to measure the likelihood of the model when it produces such predictions. For each trial, different simulations of the DEF produce different numbers of clusters, but the reported cluster parameters have fixed dimensionality for each trial. To our knowledge, modeling variable-dimensional data is not well explored in cognitive science, so we regard this as a challenge for modeling. \n\nReply to ***Is the aleatoric component definitely needed?***\n\nThe aleatoric component is definitely needed, because the rational component alone could not capture the distribution of the reported data. Please see our general response for a detailed reply. \n\nReply to ***Equation (5) may be too simple***\n\nOne could also choose other conditional distributions over K for (5); however, we want to encourage the expected number of correctly predicted K to be high, rather than to encourage e.g. the expected $\\ell-1$ or $\\ell-2$ errors to be small, such as would result from a Laplace- or Gaussian-like slack distribution.\n\n***Difference between economical and exchangeable in Figure 5(a)***\n\nThis panel only provides example predictions from the two models. The major difference is in the relative weights of the clusters. The economical model has a distortion on the reported weights, and the predicted weight values are pulled to be closer to each other (more homogeneous). \n\n***Modeling human forgetting previous examples and upweighting examples near the end***\n\nThese effects are known as recency effects in memory literature. We tested both the recency and primacy effects (up-weighting very first examples) by adding a weighting function (decaying for the primacy effect, increasing for the recency effect) to the cluster assignments. However, we did not find this to improve the model fit. We also tried to directly look for traits of primacy and recency effects in the data but did not find any significant results. \n\n***More discussions on limitations***\n\nWe will incorporate more limitations mentioned by the reviewers in a revision. In particular, we will address limitations in our experimental design, such as not asking for intermediate reports from participants. Finding a suitable way to probe for this without inducing additional cognitive effects, such as confirmation bias and forced certainty, is an important future direction. In terms of modelling, we will point out other more advanced versions of CRP, such as the repulsive CRP [\\*] that may help capture the non-overlapping clusters reported by participants. \n\n[\\*] Quinlan J J, Quintana F A, Page G L. On a class of repulsive mixture models[J]. Test, 2021\n"}]

Author rebuttal:
We thank all reviewers for their careful read of the paper and helpful feedback. We are glad that all five reviewers gave fair and comprehensive summaries, indicating that the paper is mostly well written, as also stated by four Rs (Ds9y, Xmkb, ti72 & 8bgA).  **All** Rs acknowledged that this work is interesting or novel. Three Rs (**Ds9y**, **Xmkb** & **ti72**) also praised the modelling work as clean, rigorous, non-trivial, and valuable.  

Meanwhile, most of the concerns by the Rs can be resolvable with more modelling findings (2 alternative priors for batch ICP), and richer comparisons between human behaviors and model predictions (11 figures in PDF), which we provide in this rebuttal and will incorporate in a revision. We believe these will resolve most of the weaknesses and questions of the Rs and improve the quality of the manuscript.

***How well does the model produce human data? (Ds9y, Xmkb, ti72)*** 

We initially wrote “the prediction errors on the number of clusters and cluster variance are still large” in the last paragraph as one of the limitations. We now found this to be overly pessimistic. More comprehensive comparisons between the human data and model predictions (see response PDF) indicate:
1. Our economical DEF model produces human patterns in many aspects, not only in the measures the loss function was optimized for (i.e., the number, mean, sd, and weight of clusters), but also in the moments of the whole distribution (mean, sd, and skewness) and the co-variance of different measures (e.g., how the reported cluster sd decreases with cluster number and eccentricity, but increases with cluster weight). 
2. The economical DEF performs roughly as well as the participants themselves in predicting the number of clusters $K$. In experiment 2, there are repeated trials with identical stimuli ($\mathbf{x}_T$). Using these trials, we show in Panel C of the PDF that in repeated trials participants reported a different number of clusters ($K$) just below 50% of the time. Similarly, our model could predict $K$ with accuracy around 0.5, very close to the participants themselves. The predictive power of the best model is thus limited by the stochasticity of human behavior.

***Can the data be explained by simpler models? (Ds9y, Xmkb, 8bgA)*** 

The complexity of our economical DEF as a model of our participants is not as high as it appears to be. What we assume for the participant is that they only use a *single* latent particle *at each time step* to assign the new sample to a cluster. The large number of Monte Carlo simulations is only needed for us as experimenters to fit the model, to marginalise out the latent cluster assignments in the participant’s mind and unobservable to us.

Conceptually, the CRP-GMM ICP is consistent with a simple heuristic. First, whether a new cluster should be introduced depends on a) the prior tendency to add a new cluster, described in Eqn 1; and b) how well the incoming sample is captured by the current distribution, measured by the Student’s t distribution in lines 177-178. Second, the Gaussian distribution only requires the mean (location) and variance (width), which the participants specify for each cluster.

To further validate whether the key computations in the CRP-GMM ICP are indeed necessary, we compared 11 simpler models in the paper (Appendix B.4.1). Below are examples that implement a static or less adaptive prior: 
1. The batch baseline replaces the sequential prior with a Poisson distribution and changes the structure of the ICP to a static prior;  
2. The No Counting Prior ablation replaces the counts of the clusters (core feature of CRP) with the average count; 
3. No Distribution Prior ablation removes the conjugate prior over the Gaussian distribution when evaluating the likelihood of the new sample; 

All these simpler models (see Appendix B.4.1 for the full list) led to worse fits than CRP-GMM ICP models (the exchangeable DEF and the economical DEF).

***Intuition and motivation for the aleatoric component (8bgA, ziH7)*** 

The core motivation for having the aleatoric component is to provide a well-defined likelihood function that is also cognitively plausible. Though the rational component induced by the ICP is stochastic and can in theory generate almost all behaviors with nonzero probability, we found empirically that ICP alone cannot practically predict the correct K in all trials for all subjects (i.e., likelihood is near zero). The issue lies in (as 8bgA commented) that subjects may not report exactly the clusters from their inference, due to memory and motor noises. The aleatoric component provides a highly-structured noise model to accommodate such noises (Eqns 5 & 6), which allows our model to better capture human data.  

In its technical challenge and solution, the aleatoric component resembles the classic example of drift-diffusion models (DDMs): the behavioral quantities are multi-dimensional, mixed with discrete and continuous variables; to jointly model these quantities, one can first model the discrete distributions, and then construct a conditional distribution for the continuous variables given the discrete. For DDMs, the continuous response time distribution is often conditioned on a discrete correct/wrong label. Similarly, our DEF places a conditional distribution over the discrete $K$. The rest part of the DEF is more complicated than DDMs, because 
1. the support of the predicted continuous variables (w, m, sd) depends strongly on the number of predicted clusters, and 
2. the permutation invariance of the parameters (w, m, sd) in defining the density function.

Meta-review: This paper studies how humans estimate probability distributions and presents several studies where human subject observe iid samples and are asked to estimate a mixture model distribution. The paper hypothesizes that the inconsistency between the human models and the ground truth is due to limited memory. The paper then models this process using a nonparametric Bayesian model that combines a rational component with an aleatoric component, and then evaluates this model on the data tasks.

Reviewers found the paper to contain novel contributions including a behavioral task and a new computational model. During the discussion phase, reviewers appreciated the authors rebuttal and found it clarified many questions. We encourage the authors to revise the manuscript in order to clarify many of the questions that arose during the discussion phase, including the three main questions discussed in the rebuttal below.
Decision: Accept (poster)