\documentclass[a4paper,man,natbib,floatsintext,utf8,12pt]{article}
\usepackage{geometry}
\geometry{
a4paper,
total={6in, 8in},
left=30mm,
top=30mm,
}
\renewcommand{\baselinestretch}{1.3} % spacing
\usepackage[english]{babel}
\usepackage[utf8x]{inputenc}
\usepackage{natbib}
\bibliographystyle{apalike}
\usepackage{amsmath}
\usepackage{graphicx}
\usepackage{epstopdf}
\usepackage{csquotes}
\usepackage{tabularx}
\usepackage{booktabs}
\usepackage{float}
\usepackage{amssymb}
\usepackage[dvipsnames]{xcolor}
% revision color
\newcommand{\rev}[1]{{\color{black}#1}}
% SI Figure caption
\renewcommand{\thepage}{S\arabic{page}}
\renewcommand{\thesection}{S\arabic{section}}
\renewcommand{\thetable}{S\arabic{table}}
\renewcommand{\thefigure}{S\arabic{figure}}
%linenumbers
\usepackage[section]{placeins}
\usepackage{lineno}
\linenumbers
\usepackage{hyperref}
%\usepackage{endfloat}
\hypersetup{
colorlinks,
linkcolor={red!80!black},
citecolor={blue!90!black},
urlcolor={blue!90!black}
}
\title{Supplementary Material: Entropy trade-offs in artistic design: A case study of Tamil \textit{kolam}}
\author{N.-Han Tran, Timothy Waring, Silke Atmaca, Bret A. Beheim}
\date{}
\begin{document}
\maketitle
\tableofcontents
\section{Approximation of entropy using richness and the Gini index}
\begin{figure}[h!]
\centering
\includegraphics[width=\textwidth,keepaspectratio]{SI_plots/entropy_envelop.jpeg}
\caption{High-resolution simulations showing the entropy distribution of a given richness and evenness. The black lines show the maximum entropy for a given number of variant. The differently coloured points represent the entropy distribution corresponding to the different number of variants. Equation \ref{eq:entropy_gini_n2} defines the curve for $n = 2$.}
\label{fig:entropy_envelop}
\FloatBarrier
\end{figure}
For any discrete probability distribution with $n$ possible outcomes, each $i$ of which occurs with probability $p_i$, we can calculate a number of information measures. In ecology, Shannon information entropy is a popular measure of biological diversity because it contains two different aspects: richness and evenness. The Shannon information entropy is a measure of the expected ``surprise'' or uncertainty, given in the discrete case by
\begin{equation}
H(p) = -\sum_{i=1}^{n} p_i\log(p_i)
\label{eq:entropy}
\end{equation}
for each outcome $i$. In economics, the Gini index \citep{Zoli1999, Ravallion2014} is used to describe relative inequality in the probability distribution, calculated as the mean absolute difference between all pair of variants,
\begin{equation}
g(n) = \frac{\sum_{i=1}^{n}\sum_{j=1}^{n} |p_i-p_j|}{2(n-1)}
\label{eq:gini}
\end{equation}
where $n$ is the richness (the number of unique variants) and $p$ the frequency of specific variants.
Entropy and the Gini index capture variation in the relative abundance of each outcome following the focal probability distribution. That is, when one outcome $p_j \to 1$ and all $p_{_j} \to 0$, the Shannon entropy goes to its lower asymptote of 0, and the Gini index goes to its upper asymptote of 1. Likewise, when $p_j = \frac{1}{n}$ for all $j$, the entropy is maximized at \rev{$\log(n)$} and the Gini index is minimized at 0.
In the special case of $n = 2$ with two variant frequencies $p$ and $q = 1 - p$, such that $p > q$, the Gini index simplifies to $g_2 = p - q$. Using the fact that $p + q = 1$, we can rewrite each as: $p = \frac{1-g_2}{2}, q = \frac{1+g_2}{2}$. Thus, there is an exact relationship between $H_2$ and $g_2$, such that entropy is maximized when the Gini index is minimized, and vice versa.
\begin{equation}
H_2 = \left(\frac{1+g}{2}\right) \log\left(\frac{2}{1+g}\right) + \left(\frac{1-g}{2}\right) \log\left(\frac{2}{1-g}\right)
\label{eq:entropy_gini_n2}
\end{equation}.
Equation \ref{eq:entropy_gini_n2} defines the curve for $n = 2$ in Figure \ref{fig:entropy_envelop}.
However, the relationship between a Gini index and an entropy is indefinite if $n > 2$ because multiple distributions with the same $n$ with different entropy could take the same Gini index value. Figure \ref{fig:entropy_envelop} illustrates the relationship between the entropy and Gini index calculated for 100,000 simulated probability distributions.
To understand the relationship between entropy, the Gini index and the richness further, the location of the minimum and maximum entropy within this wing-shaped ``envelope'' are important. \rev{Given any particular value of $g$, and number of variants $n$, we can describe the range of possible distributions between a maximum and minimum entropy. Although analytic solutions for the minimum or maximum entropy are elusive for $n > 3$, numerical solutions are readily available using simulation and non-linear optimization algorithms. In the supplementary codebase, we use the \texttt{Rsolnp} package to run the non-linear optimization algorithm and find maximum-entropy solutions for the cases in Figure \ref{fig:entropy_envelop}.}
The lower ``tips'' of each distribution in Figure \ref{fig:entropy_envelop} represent the minimum-entropy limits at which the least-common non-zero variant tends to a zero frequency. \rev{In the case of $n = 3$, the} minimum-entropy boundary can be approximated by noting that the \rev{minimum entropy (i.e., the ``tip'')} occurs when $g = \frac{1}{2}$. At a closer look at Figure \ref{fig:entropy_envelop}, \rev{the minima of the entropy distribution lie at $\{\frac{0}{n-1}, \frac{1}{n-1},\dots,\frac{n-1}{n-1}\}$, following a limiting boundary described by the equation:}
\begin{equation}
H_{min} = \log\left(n-(n - 1)g\right).
\label{eq:entropy_min}
\end{equation} \rev{More precisely,} the minimum entropy boundaries can be defined by generalizing the entropy equation for n = 2 in equation \ref{eq:entropy_gini_n2} to: $p = \frac{y - a}{b - a}, q = \frac{b - y}{b-a}$, where $v = 1 - g$ in the special case of $n = 2, a = 0$ and $b = 2$. \rev{Empirical exploration of the entropy envelope by the simulation of 100,000 gesture distributions indicates that no distribution can occupy a lower entropy than defined by this boundary.} This can be used to decompose entropy using $g$ and $n$ because if $g$ is known, $H_{n, g}$ varies between boundaries defined by the equation:
\begin{equation}
\mathrm{exp}(\widehat{H}) \approx n - (n-1) v^{1 + b} = n - (n-1) v^{1 + \frac{2}{2 + a} + \frac{a}{a + n}},
\label{eq:entropy_derivation}
\end{equation} where evenness is $v = 1 - g$\rev{. If $b = 0$, this equation gives the lower theoretical boundary on entropy for any given $g$ and $n$ described by equation \ref{eq:entropy_min} and the ``tips'' in Figure \ref{fig:entropy_envelop}. The numerically-calculated maximum entropy values are described by the equation $b = \frac{2}{2 + a} + \frac{a}{a + n}$ with $a \approx \exp(0.5139)$. Thus, if we assume that a distribution tends towards maximum entropy \citep{Frank2009}, we can calculate a distributionâ€™s entropy knowing only its richness $n$ and the evenness $v = 1 - g$.}
\rev{As evenness and richness are now related to entropy by a mathematical identity (and so have no single causal direction to their relationship), we can aspire to understand what actually could explain why the observed \emph{kolam} patterns follow the entropy isoclines.}
\section{\emph{Kolam} Data}
\subsubsection{Lexicon of Gestures}
One specific class of \emph{kolam} are those with loop patterns, called square loop \emph{kolam} drawings (i.e., the \emph{ner pulli nelevu} or \emph{sikku} \emph{kolam} family) \citep{Waring}. These \emph{kolam} drawings are composed of an initial grid of dots (\emph{pulli}) that reflect the canvas size. Lines consisting of multiple gestures are sequentially drawn around the dots to form loops. We only focus on these square loop \emph{kolam} drawings because the patterns can be mapped onto a small identifiable set of gestures which is suitable for analyses.
The geometry of the \emph{kolam} can be divided in two fundamental geometric spaces and a transitional geometric space with their specific corresponding positions, orientations and gestures. All gestures are always located in relation to the neighboring dots, called \textit{pulli}. For a detailed description of the sequential encoding of \emph{kolam} drawings with the gestural lexicon, please consult \citet{Waring}.
\begin{figure}[h!]
\centering
\includegraphics[width=\textwidth,height=0.7\textheight,keepaspectratio]{SI_plots/Lexicon.jpeg}
\label{fig:lexicon2}
\caption{The Lexicon of \emph{Kolam} Gestures. The Figure illustrates the gestures and the corresponding code to encode \emph{kolam} drawings. Taken and adapted with permission from \citet{Waring}}
\end{figure}
\FloatBarrier
\subsection{Database}
\subsubsection{Survey Information}
Information on individuals' \emph{kolam} drawing abilities and behaviour were gathered as well as demographic information. Demographic information entailed data such as GPS data of the current residency, marriage status, native places and measures of socio-economic status. To investigate individual's \emph{kolam} practice, survey questions included information on individual's frequency of drawing \emph{kolam} drawings in front of their doorstep or in their practice book and the age of initial learning.
\subsubsection{GPS}
The geographical position of each individual's current residency was measured. A GPS tracker of type Garmin GPSmap 60 CSx was used. The interviewer only asked for the name of the native place and during data pre-processing the name of the native place and if needed other demographic information was then used to manually map the name of the native place to a GPS position. All GPS positions of the current residency and the native place were recorded allowing distances between points to be calculated by applying the distance formula to the x-y-z coordinates of the two points.
\subsubsection{\emph{Kolam} Drawings}
A corpus of \emph{kolam} drawings was compiled by soliciting individual's to draw \emph{kolam} drawings as part of the survey. Each individual was asked to draw a minimum of a total of 20 \emph{kolam} drawings.
\section{Statistical Analyses}
\subsection{Random Intercept Models}
\subsubsection{Statistical Model}
\begin{equation}
\label{eq1}
\begin{split}
\text{Density} \sim \text{Log-Normal}(\mu_i, \sigma) \\
\mu_i = \alpha_{density} + \alpha_j + \beta_{caste} + \beta_{age} \times \text{age} + \beta_{\text{practice duration}} \times \text{practice duration} \\
\beta_{caste} = \sigma_{caste} \times z_{caste} \\
\alpha_{j} = \sigma_{artist} \times z_{artist} \\
\alpha_{density} \sim \text{Normal}(1, 2) \\
\sigma_{caste} \sim \text{Normal}(0.5, 1) \\
\sigma_{artist} \sim \text{Normal}(0, 0.5) \\
\sigma \sim \text{Normal}(0.5, 0.5) \\ \beta_{caste} \sim \text{Normal}(0, 1) \\
\beta_{age} \sim \text{Normal}(0, 1) \\
z_{caste} \sim \text{Normal}(0, 1) \\
z_{artist} \sim \text{Normal}(0, 1)
\end{split}
\end{equation}
\begin{equation}
\label{eq2}
\begin{split}
\text{Evenness} \sim \text{Truncated Normal}(\mu_i, \sigma)[0, 1] \\
\mu_i = \alpha_{evenness} + \alpha_j + \beta_{caste} + \beta_{age} \times \text{age} + \beta_{\text{practice duration}} \times \text{practice duration} \\
\beta_{caste} = \sigma_{caste} \times z_{caste} \\
\alpha_{j} = \sigma_{artist} \times z_{artist} \\
\alpha_{evenness} \sim \text{Normal}(1, 2) \\
\sigma_{caste} \sim \text{Normal}(0.5, 1) \\
\sigma_{artist} \sim \text{Normal}(0, 0.5) \\
\sigma \sim \text{Normal}(0.5, 1) \\ \beta_{caste} \sim \text{Normal}(0, 1) \\
\beta_{age} \sim \text{Normal}(0, 1) \\
z_{caste} \sim \text{Normal}(0, 1) \\
z_{artist} \sim \text{Normal}(0, 1)
\end{split}
\end{equation}
\begin{equation}
\label{eq3}
\begin{split}
\text{Richness} \sim \text{Poisson}(\lambda_i) \\
\log(\lambda_i) = \alpha_{richness} + \alpha_j + \beta_{caste} + \beta_{age} \times \text{age} + \beta_{\text{practice duration}} \times \text{practice duration} \\
\beta_{caste} = \sigma_{caste} \times z_{caste} \\
\alpha_{j} = \sigma_{artist} \times z_{artist} \\
\alpha_{richness} \sim \text{Normal}(1, 2) \\
\sigma_{caste} \sim \text{Normal}(0.5, 1) \\
\sigma_{artist} \sim \text{Normal}(0, 0.5) \\
\beta_{caste} \sim \text{Normal}(0, 1) \\
\beta_{age} \sim \text{Normal}(0, 1) \\
z_{caste} \sim \text{Normal}(0, 1) \\
z_{artist} \sim \text{Normal}(0, 1)
\end{split}
\end{equation}
\newpage
\begin{equation}
\label{eq4}
\begin{split}
\text{Canvas Size}_{i} \sim \text{NegBinom}(\mu_i, \phi) \\
\log(\mu_i) = \alpha_{size} + \alpha_{j} + \beta_{caste} + \beta_{age} \times \text{age} + \beta_{\text{practice duration}} \times \text{practice duration} \\
\beta_{caste} = \sigma_{caste} \times z_{caste} \\
\alpha_{j} = \sigma_{artist} \times z_{artist} \\
\alpha_{size} \sim \text{Normal}(1, 2) \\
\sigma_{caste} \sim \text{Normal}(0.5, 1) \\
\sigma_{artist} \sim \text{Normal}(0, 0.5) \\
\phi \sim \text{Normal(1.5, 3)} \\
\beta_{caste} \sim \text{Normal}(0, 1) \\
\beta_{age} \sim \text{Normal}(0, 1) \\
z_{caste} \sim \text{Normal}(0, 1) \\
z_{artist} \sim \text{Normal}(0, 1)
\end{split}
\end{equation}
\begin{equation}
\label{eq5}
\begin{split}
\text{Entropy}_{i} \sim \text{Truncated Normal}(\mu_i, \sigma)[0, 1] \\
\mu_i = \alpha_{entropy} + \alpha_{j} + \beta_{caste} + \beta_{age} \times \text{age} + \beta_{\text{practice duration}} \times \text{practice duration} \\
\beta_{caste} = \sigma_{caste} \times z_{caste} \\
\alpha_{j} = \sigma_{artist} \times z_{artist} \\
\alpha_{entropy} \sim \text{Normal}(1, 2) \\
\sigma_{caste} \sim \text{Normal}(0.5, 1) \\
\sigma_{artist} \sim \text{Normal}(0, 0.5) \\
\phi \sim \text{Normal(1.5, 3)} \\
\beta_{caste} \sim \text{Normal}(0, 1) \\
\beta_{age} \sim \text{Normal}(0, 1) \\
z_{caste} \sim \text{Normal}(0, 1) \\
z_{artist} \sim \text{Normal}(0, 1)
\end{split}
\end{equation}
\FloatBarrier
\rev{
\subsubsection{Visual MCMC Diagnostics}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/trace_density_intercept.png}
\caption{Traceplot for the random intercept model on density showing mixing across chains and convergence.}
\label{fig:trace_density}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/trace_gini_intercept.png}
\caption{Traceplot for the random intercept model on evenness showing mixing across chains and convergence.}
\label{fig:trace_evenness}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.70\textwidth,keepaspectratio]{SI_plots/trace_unique_intercept.png}
\caption{Traceplot for the random intercept model on richness showing mixing across chains and convergence.}
\label{fig:trace_richness}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/trace_size_intercept.png}
\caption{Traceplot for the random intercept model on canvas size showing mixing across chains and convergence.}
\label{fig:trace_size}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/trace_entropy_intercept.png}
\caption{Traceplot for the random intercept model on entropy showing mixing across chains and convergence.}
\label{fig:trace_entropy}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/pairs_density_intercept.png}
\caption{Pairs plot for the random intercept model on density showing correlation among parameters and no sampling problems.}
\label{fig:pairs_density}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/pairs_gini_intercept.png}
\caption{Pairs plot for the random intercept model on gini showing correlation among parameters and no sampling problems.}
\label{fig:pairs_evenness}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.70\textwidth,keepaspectratio]{SI_plots/pairs_unique_intercept.png}
\caption{Pairs plot for the random intercept model on richness showing correlation among parameters and no sampling problems.}
\label{fig:pairs_richness}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/pairs_size_intercept.png}
\caption{Pairs plot for the random intercept model on canvas size showing correlation among parameters and no sampling problems.}
\label{fig:pairs_size}
\end{figure}
\begin{figure}[!ht]
\centering
\includegraphics[width=0.75\textwidth,keepaspectratio]{SI_plots/pairs_entropy_intercept.png}
\caption{Pairs plot for the random intercept model on entropy showing correlation among parameters and no sampling problems.}
\label{fig:pairs_entropy}
\end{figure}
}
\FloatBarrier
\subsubsection{Intraclass Correlation (ICC)}
The intraclass-correlation coefficient (ICC) can be calculated for Gaussian models to determine the variance explained by random and fixed effects \citep[p.258]{gelman_hill_2006}. Since our five models are non-Gaussian, we approximated the ICC using variance decomposition based on the posterior predictive distribution. We first drew from the posterior predictive distribution not conditioned on our fixed (i.e., age and duration of practice) and random effect (i.e., caste and individual variation) terms and then drew from the posterior predictive distribution conditioned on all fixed and random effects. Subsequently, we calculated the variances for both draws. The pseudo-ICC is then the ratio between these two variances. Occasionally, the variance ratio can be negative due to very large variance of the posterior predictive distributions.
\begin{figure}[!htbp]
\centering
\includegraphics[width=0.7\textwidth,keepaspectratio]{SI_plots/ICC.jpeg}
\caption{Intraclass Correlation Coefficients (ICC) for individual random-effect variances for the four outcome variables. The variance decomposition is based on the posterior predictive distribution, which is the correct way for Bayesian non-Gaussian models.}
\label{fig:ICC}
\end{figure}
\FloatBarrier
\subsubsection{Illustration of Random Effect Estimates for Caste}
\FloatBarrier
\begin{figure}[!htbp]
\centering
\includegraphics[width=\textwidth,height=\textheight,keepaspectratio]{SI_plots/caste_offset_sort.eps}
\caption{Caste random effect offsets for each outcome variable. The red line reflects zero offset. Each violin probability density plot displays the variation within each caste on the corresponding outcome variable (i.e., richness, canvas size, density, evenness and entropy). The posterior mean offset is illustrated in blue. The range of the violin plots reflect the $90\%$ HPDI.}
\label{fig:caste_offset}
\end{figure}
\FloatBarrier
\subsection{Gaussian Process Models}
To investigate whether the observed variation in structural and information-theoretic properties of \emph{kolam} drawings is structured by the residential or native place proximity between individuals, we used a Gaussian process (GP) model to estimate a function for the covariance between pairs of individuals at different spatial distances of their residency as well as their native place. We \rev{fit} the GP model on the five outcome variables entropy, canvas size, gesture density, richness and evenness. To foreshadow the conclusions drawn from the GP model, results from the GP model are in line with our main results from the random-intercept models. The GP model further indicates no major differentiation between artists originally from the community and those who emigrated from other parts of India, and no spatial association between artist's \emph{kolam} drawings.
\subsubsection{Statistical Model}
We estimated unique intercepts for each individual with a varying effects approach to continuous categories using a Gaussian process (GP) model. A GP can be seen as a distribution of nonlinear functions. Placing a GP prior over the covariance, allows us to estimate a function for the covariance between pairs of individuals at different variable distances \citep{McElreath2016}. The individual-level covariance matrix $\alpha_{i}$ was modeled using an exponentiated quadratic kernel. This covariance function implies that the covariance between any two individuals $j$ and $k$ declines exponentially with the squared distance between them. No variation between individuals would correspond to a covariance of zero or close to zero. The parameter $\rho^2$, also referred to as length scale, determines the rate of decline. If $\rho^2$ is large, then the covariance decreases slowly with squared distance, while if $\rho^2$ is small it decreases rapidly with squared distance. This length-scale $\rho^2$ prior is constrained to be between zero and one because the distances were normalized to be between zero and one. $\eta^2$ is the maximum covariance between any two individuals. Automatic relevance determination \citep{Neal1996} was performed where multiple predictors with corresponding length-scale parameters for each dimension were merged into the GP.
Statistical models on the five structural and information-theoretic properties of \emph{kolam} drawings were implemented. The predictor variables to investigate the partition of variation across individuals and the population was the same in all models. Age, duration of practice, caste, and residential and native place proximity between individuals were used as predictor variables to explain individual variation on the structural and information-theoretic properties. The residential and native place distances as well as age and practice duration were normalized, such that the minimum value was zero and the maximum value was one. The pairwise Euclidean distances between each pair of individuals along each predictor variable dimension was then computed for age, duration of practice and residential and native place proximity between individuals. GPS coordinates were used to compute residential and native place distance matrices between individuals. Distances correspond to Euclidean distance between individuals for the specific variable dimension. For example, the residence or native place distance matrix are spatial distances. Caste was modeled as a hierarchical non-centred categorical variable with 19 categories.
\begin{equation*}
\label{gp_eq1}
\begin{split}
\text{Density} \sim \text{Log-Normal}(\mu_i, \sigma) \\
\mu_i = \alpha_{density} + \alpha_j + \beta_{caste}\\
\alpha_j \sim \text{MVNormal}(0, K(x)) \\
K_{jk} = \eta^2 \text{exp}\left[-\left(\frac{\text{Residence}^2}{2 \times \rho_R^2} + \frac{\text{Native Place}^2}{2 \times \rho_N^2} + \frac{\text{Age}^2}{2 \times \rho_A^2} + \frac{\text{Duration}^2}{2 \times \rho_D^2}\right)\right] + \delta_{jk} \times 0.001 \\
\beta_{caste} = \sigma_{caste} \times z \\
\eta^2 \sim \text{Normal}(5, 2) \\
\rho_R^2 \sim \text{Beta}(1, 2) \\
\rho_N^2 \sim \text{Beta}(1, 2) \\
\rho_A^2 \sim \text{Beta}(1, 2) \\
\rho_D^2 \sim \text{Beta}(1, 2) \\
\sigma \sim \text{Normal}(0.5, 0.5) \\
\alpha_{density} \sim \text{Normal}(0.5, 1) \\
\sigma_{caste} \sim \text{Normal}(0.5, 0.5) \\
\end{split}
\end{equation*}
\begin{equation*}
\label{gp_eq2}
\begin{split}
\text{Evenness} \sim \text{Truncated Normal}(\mu_i, \sigma)[0, 1] \\
\mu_i = \alpha_{evenness} + \alpha_j + \beta_{caste} \\
\alpha_j \sim \text{MVNormal}(0, K(x)) \\
K_{jk} = \eta^2 \text{exp}\left[-\left(\frac{\text{Residence}^2}{2 \times \rho_R^2} + \frac{\text{Native Place}^2}{2 \times \rho_N^2} + \frac{\text{Age}^2}{2 \times \rho_A^2} + \frac{\text{Duration}^2}{2 \times \rho_D^2}\right)\right] + \delta_{jk} \times 0.001 \\
\beta_{caste} = \sigma_{caste} \times z \\
\eta^2 \sim \text{Normal}(5, 2) \\
\rho_R^2 \sim \text{Beta}(1, 2) \\
\rho_N^2 \sim \text{Beta}(1, 2) \\
\rho_A^2 \sim \text{Beta}(1, 2) \\
\rho_D^2 \sim \text{Beta}(1, 2) \\
\sigma \sim \text{Normal}(0.5, 0.5) \\
\alpha_{evenness} \sim \text{Normal}(0.5, 1) \\
\sigma_{caste} \sim \text{Normal}(0.5, 0.5) \\
\end{split}
\end{equation*}
\begin{equation*}
\label{gp_eq3}
\begin{split}
\text{Richness} \sim \text{Poisson}(\lambda_i) \\
\log(\lambda_i) = \alpha_{richness} + \alpha_j + \beta_{caste} \\
\alpha_j \sim \text{MVNormal}(0, K(x)) \\
K_{jk} = \eta^2 \text{exp}\left[-\left(\frac{\text{Residence}^2}{2 \times \rho_R^2} + \frac{\text{Native Place}^2}{2 \times \rho_N^2} + \frac{\text{Age}^2}{2 \times \rho_A^2} + \frac{\text{Duration}^2}{2 \times \rho_D^2}\right)\right] + \delta_{jk} \times 0.001 \\
\beta_{caste} = \sigma_{caste} \times z \\
\eta^2 \sim \text{Normal}(5, 2) \\
\rho_R^2 \sim \text{Beta}(1, 2) \\
\rho_N^2 \sim \text{Beta}(1, 2) \\
\rho_A^2 \sim \text{Beta}(1, 2) \\
\rho_D^2 \sim \text{Beta}(1, 2) \\
\sigma \sim \text{Normal}(0.5, 0.5) \\
\alpha_{richness} \sim \text{Normal}(0.5, 1) \\
\sigma_{caste} \sim \text{Normal}(0.5, 0.5) \\
\end{split}
\end{equation*}
\newpage
\begin{equation*}
\label{gp_eq4}
\begin{split}
\text{Canvas Size}_{i} \sim \text{NegBinom}(\mu_i, \phi) \\
\log(\mu_i) = \alpha_{size} + \alpha_{j} + \beta_{caste} \\
\alpha_j \sim \text{MVNormal}(0, K(x)) \\
K_{jk} = \eta^2 \text{exp}\left[-\left(\frac{\text{Residence}^2}{2 \times \rho_R^2} + \frac{\text{Native Place}^2}{2 \times \rho_N^2} + \frac{\text{Age}^2}{2 \times \rho_A^2} + \frac{\text{Duration}^2}{2 \times \rho_D^2}\right)\right] + \delta_{jk} \times 0.001 \\
\beta_{caste} = \sigma_{caste} \times z \\
\eta^2 \sim \text{Normal}(5, 2) \\
\rho_R^2 \sim \text{Beta}(1, 2) \\
\rho_N^2 \sim \text{Beta}(1, 2) \\
\rho_A^2 \sim \text{Beta}(1, 2) \\
\rho_D^2 \sim \text{Beta}(1, 2) \\
\sigma \sim \text{Normal}(0.5, 0.5) \\
\alpha_{size} \sim \text{Normal}(0.5, 1) \\
\sigma_{caste} \sim \text{Normal}(0.5, 0.5) \\
\end{split}
\end{equation*}
\begin{equation*}
\label{gp_eq5}
\begin{split}
\text{Entropy} \sim \text{Truncated Normal}(\mu_i, \sigma)[0, ] \\
\mu_i = \alpha_{entropy} + \alpha_j + \beta_{caste} \\
\alpha_j \sim \text{MVNormal}(0, K(x)) \\
K_{jk} = \eta^2 \text{exp}\left[-\left(\frac{\text{Residence}^2}{2 \times \rho_R^2} + \frac{\text{Native Place}^2}{2 \times \rho_N^2} + \frac{\text{Age}^2}{2 \times \rho_A^2} + \frac{\text{Duration}^2}{2 \times \rho_D^2}\right)\right] + \delta_{jk} \times 0.001 \\
\beta_{caste} = \sigma_{caste} \times z \\
\eta^2 \sim \text{Normal}(5, 2) \\
\rho_R^2 \sim \text{Beta}(1, 2) \\
\rho_N^2 \sim \text{Beta}(1, 2) \\
\rho_A^2 \sim \text{Beta}(1, 2) \\
\rho_D^2 \sim \text{Beta}(1, 2) \\
\sigma \sim \text{Normal}(0.5, 0.5) \\
\alpha_{entropy} \sim \text{Normal}(0.5, 1) \\
\sigma_{caste} \sim \text{Normal}(0.5, 0.5) \\
\end{split}
\end{equation*}
\subsubsection{Estimation of Variation}
\FloatBarrier
The five statistical models were implemented in the probabilistic programming language Stan 2.18 \citep{Stan}, using 6000 samples from four chains. Analyses were performed in R \citep{RCoreTeam2019}. Data and analyses can be found here: \url{http://github.com/nhtran93/kolam_signaling}. All R-hat values were less than 1.01, and visual inspection of trace plots and rank histograms indicated convergence of all models.
Between-individual variation in entropy ($\eta^2$ = 0.00, $90\%$ CI [0.00, 0.00]), density ($\eta^2$ = 0.01, $90\%$ CI [0.01, 0.01]), the evenness ($\eta^2$ = 0, $90\%$ CI [0, 0]), and in the richness ($\eta^2$ = 0.00, 90\% CI [0, 0.01]) were estimated with high certainty to be very small and close to zero (see left panel in Figure \ref{fig:coefficients}). The between-artist variability is most pronounced in canvas size (($\eta^2$ = 0.03, $90\%$ CI [0.02, 0.03]), while \emph{kolam} drawings show only small distinct variation between artists in the other structural and information-theoretic properties. Our results further show no evidence that variation in structural and information-theoretic properties of \emph{kolam} drawings covary with the spatial structure, age or practice duration (see Figure \ref{fig:coefficients}).
Prior-posterior plots of all five models show that the priors updated for all parameters except the length scale parameters $\rho$ because there is barely any information to explain individual variation with no individual variation present. We detected very small effects of caste membership on the entropy, density, evenness, and richness, with varying-effect deviations estimated near zero with high certainty (entropy: $\sigma_{caste}$ = 0.03, $90\%$ CI [0.01, 0.05]; density: $\sigma_{caste}$ = 0.04, $90\%$ CI [0.01, 0.07]; evenness: $\sigma_{caste}$ = 0.03, $90\%$ CI [0.01, 0.04]; and richness: $\sigma_{caste}$ = 0.01, $90\%$ CI [0.00, 0.03] respectively). We detected more pronounced effects of caste membership on canvas size (canvas size: $\sigma_{caste}$ = 0.09, $90\%$ CI [0.05, 0.14]) as illustrated in Figure \ref{fig:coefficients}.
\begin{figure} [h!]
\centering
\includegraphics[width=\textwidth,height=0.9\textheight,keepaspectratio]{SI_plots/coefficients_GP.eps}
\caption{Prior-Posterior Coefficient Plots of Individual variation and Caste Variation. All panels have the same y-axis indicating the five models. The left panel (eta squared) illustrates the estimated individual variation (dark blue) in comparison to the prior (light blue) for each model. The right panel illustrates the estimated population-level standard deviation for the effect of caste (dark blue) in comparison to the prior (light blue) for each model. The 90\% Highest Posterior Density Interval (HPDI) was computed for each posterior; however, the interval is very narrow.}
\label{fig:coefficients}
\end{figure}
\bibliography{bib}
\end{document}