Intermittency Analysis of Toy Monte Carlo Events

Sheetal Sharma and Ramni Gupta$\star$

Department of Physics, University of Jammu, India

*<EMAIL_ADDRESS>
50th International Symposium on Multiparticle Dynamics
(ISMD2021)
12-16 July 2021 10.21468/SciPostPhysProc.?

## Abstract

Event-by-event intermittency analysis of Toy Monte Carlo events is performed
in the scenario of high multiplicity events as is the case at recent colliders
RHIC and LHC for AA collisions. A power law behaviour of Normalized Factorial
Moments (NFM), $F_{q}$ as function of number of bins ($M$) known as
intermittency, is a signature of self-similar fluctuations. Dependence of NFM
on the detector efficiencies and on the presence of fluctuations have been
studied. Results presented here provide a baseline to the experimental results
and clarity on the application of efficiency corrections to the experimental
data.

##  1 Introduction

Localized fluctuations in the charged particle production at LHC are proposed
to be studied to characterize the multiparticle production and the quark-
hadron phase transition [1]. QCD predicts large dynamical fluctuations of
various measureables as one of the signatures of critical point, quark-hadron
and hadron-quark phase transition. A study of spatial patterns of the charged
particles in the phase space using normalized factorial moments (NFM) is one
of the techniques to characterize phase transition and the multiparticle
production mechanism [2, 3, 4]. Normalized factorial moments ($F_{q}$) of bin
multiplicities as function of varying bin size resolution are proposed to be
studied for q $\geq$ 2 [1]. For dynamical fluctuations $F_{q}>1$ and is
observed to show power-law behaviour with increasing M for self similar
fluctuations and this phenomenon is known as intermittency. Here intermittency
analysis is performed for the Toy Monte Carlo (ToyModel) events as baseline
study.

##  2 Method of Analysis

A sample of 250K high multiplicity Toy Monte Carlo events are generated with
two parameters for the tracks corresponding to pseudorapidity $(\eta)$ and
azimuthal angle $(\phi)$ such that $\mid\eta\mid\leq 0.8$ and $0\leq\phi\leq
6.28$. Intermittency analysis (as in [5]) is performed in two dimensional
$(\eta,\phi)$ phase space partitioned in $M\times M$ bins, with M = 4 to 82.
The $q^{th}$ order normalized factorial moment ($F_{q}$) is defined as

$F_{q}(M)=\frac{\displaystyle{\frac{1}{N}}\displaystyle\sum_{e=1}^{N}\frac{1}{M}\sum_{i=1}^{M}f_{q}(n_{ie})}{\left(\displaystyle{\frac{1}{N}}\displaystyle\sum_{e=1}^{N}\frac{1}{M}\sum_{i=1}^{M}f_{1}(n_{ie})\right)^{q}}$
(1)

where $f_{q}(n_{ie})=\Pi_{j=0}^{q-1}(n_{ie}-j)$, $n_{ie}$ is the bin
multiplicity in the $i^{th}$ bin of $e^{th}$ event. $q\geq 2$ is order of the
moment and takes positive integer values. $F_{q}(M)$ shows power law
dependence on M as $F_{q}(M)\propto M^{\phi_{q}}$ with $\phi_{q}>0$ in case
there are fluctuations in the bin multiplicities [5]. This scaling behaviour
is referred to as  intermittency and $\phi_{q}$ as intermittency index. With
second order NFM (q = 2), the sensitivity of this analysis methodology to
gauge bin-to-bin fluctuations and the resilience to detector inefficiencies
has been studied in the present work.

##  3 Observations

Normalized factorial moments for q = 2 are determined for Toy Monte Carlo
events (ToyModel) using Eq.1. It is observed that for all M, $F_{2}(M)>1$
(Fig.1, black filled circles). Also $F_{2}$ values are independent of number
of bins (M). Toy Monte Carlo events do not show any scaling behaviour and
hence no intermittency.

For sensitivity check of the analysis methodology a modified sample of events
is created from the ToyModel events, using two different approaches. In the
first method five percent tracks are added randomly in some phase space bins
and an equal number of tracks are removed from rest of the region. In the
second method, in a similar fashion five percent tracks are added in some
phase space bins but no tracks are removed. For both samples so obtained,
intermittency analysis is performed and it is observed that $F_{2}>1$. However
$F_{2}$ is observed to depend on M (Fig.1). At higher M region $F_{2}$ shows
has linear growth with M. This establishes that intermittency analysis
methodology is sensitive to the particle density fluctuations.

Calculations for the observables are affected by the detector effects and
hence do not give true value. If $\epsilon_{i}$ defines the detector
efficiency in the $i^{th}$ bin then corrected NFM is taken as

$F_{q}(M)=\frac{\displaystyle{\frac{1}{N}}\displaystyle\sum_{e=1}^{N}\frac{1}{M}\sum_{i=1}^{M}\frac{f_{q}(n_{ie})}{\epsilon_{i}^{q}}}{\left(\displaystyle{\frac{1}{N}}\displaystyle\sum_{e=1}^{N}\frac{1}{M}\sum_{i=1}^{M}\frac{f_{1}(n_{ie})}{\epsilon_{i}}\right)^{q}}$
(2)

From the ToyModel events, which may be called ToyModel(true), two samples of
events equivalent to what is measured by the detectors after applying
reconstruction routine are obtained. First sample is created by randomly
removing 20% of tracks from the acceptance region of each event. The sample so
obtained, say ToyModel(U), is 80% of the ToyModel(true) and has uniform
efficiency across the acceptance region. For events with non-binomial type
efficiencies, 20% particles are removed from some specific phase space regions
of each event. This sample of events, say ToyModel(NU), is also 80% of the
original events but with different efficiencies across the acceptance region.

Figure 1: $F_{2}$ vs $M^{2}$, depicting sensitivity of analysis technique to
gauge bin-to-bin fluctuations in the ToyModel events.

Normalized factorial moments are determined for the three samples of events
using Eq.1 and corrected NFM are determined for ToyModel(U) and ToyModel(NU)
using Eq.2. It is observed that $F^{(true)}_{2}(M)\approx
F^{(U)}_{2}(M)\approx F^{(Ucorr)}_{2}(M)$ as is shown in Fig.2(a). However as
in Fig.2(b) $F^{(true)}_{2}(M)\neq F^{(NU)}_{2}(M)$ (black open circle and
blue solid square markers), that is any change in the true track values,
introduced differently in different phase space regions, the NFM (Eq.1) do not
give true value. Whereas $F^{(true)}_{2}(M)\approx F^{(NUcorr)}_{2}(M)$
(Fig.2(b)) implying that with corrected NFM calculated using Eq.2 the true NFM
are reproduced. Thus the NFM are robust against binomial detector efficiencies
but for non-binomial detector efficiencies, to obtain true NFM, formula (Eq.2)
with bin efficiency correction values must be used.

Figure 2: $F_{2}$ vs $M^{2}$ plot in case of (a)binomial type efficiencies and
(b)non-binomial type efficiencies.

##  4 Conclusions

Intermittency analysis is performed for high multiplicity Toy Monte Carlo
events. Analysis technique is observed to be suitable to look for dynamical
fluctuations in the multiplicity distributions. NFM as defined in Eq.1 are
observed to be robust against the binomial detector efficiencies. However NFM
should be corrected for detector effects if efficiencies are non-binomial/non-
Gaussian in the acceptance region before any conclusions be drawn.

## Acknowledgements

Author(s) are thankful to Igor Altsybeev, Mesut Arslandok and Tapan Nayak for
the discussions and suggesstions that helped in the completion of this study.

### Author contributions

R.G. designed the model and the computational framework. Both R.G. and S.S.
performed the analysis and worked to bring out this manuscript.

### Funding information

The author(s) received no financial support for the research, authorship
and/or publication of this work.

## References

  * [1] R. C. Hwa and C. B. Yang, Local Multiplicity Fluctuations as a Signature of Critical Hadronization at LHC, Phys. Rev. C 85, 044914 (2012), 10.1103/PhysRevC.85.044914.
  * [2] A. Bialas and R. Peschanski, Moments of rapidity distributions as a measure of short-range ﬂuctuations in high-energy collisions, Nucl. Phys. B 273, 703 (1986), 10.1016/0550-3213(86)90386-X.
  * [3] R. C. Hwa and C. B. Yang, Observable Properties of Quark-Hadron Phase Transition at the Large Hadron Collider, Acta Phys. Polon. B 48, 23 (2017), 10.5506/APhysPolB.48.23.
  * [4] W. Kittel and E. A. De Wolf, Soft Multihadron Dynamics, World Scientiﬁc, Singapore, 2005.
  * [5] R. Sharma and R. Gupta, Scaling Properties of Multiplicity Fluctuations in the AMPT Model, Adv. High Energy Phys. 2018, 6283801 (2018), 10.1155/2018/6283801.


# Pointer-Generator Networks for Low-Resource Machine Translation: Don’t Copy
That!

Niyati Bafna, Philipp Koehn, and David Yarowsky
Johns Hopkins University, Center for Language and Speech Processing
<EMAIL_ADDRESS>
###### Abstract

While Transformer-based neural machine translation (NMT) is very effective in
high-resource settings, many languages lack the necessary large parallel
corpora to benefit from it. In the context of low-resource (LR) MT between two
closely-related languages, a natural intuition is to seek benefits from
structural “shortcuts”, such as copying subwords from the source to the
target, given that such language pairs often share a considerable number of
identical words, cognates, and borrowings. We test Pointer-Generator Networks
for this purpose for six language pairs over a variety of resource ranges, and
find weak improvements for most settings. However, analysis shows that the
model does not show greater improvements for closely-related vs. more distant
language pairs, or for lower resource ranges, and that the models do not
exhibit the expected usage of the mechanism for shared subwords. Our
discussion of the reasons for this behaviour highlights several general
challenges for LR NMT, such as modern tokenization strategies, noisy real-
world conditions, and linguistic complexities. We call for better scrutiny of
linguistically motivated improvements to NMT given the blackbox nature of
Transformer models, as well as for a focus on the above problems in the field.

##  1 Introduction and Motivation

While state-of-the-art Transformer models (Vaswani et al., 2017) for NMT work
well for high-resource language pairs, their performance degrades in low-
resource situations (Kim et al., 2020; Haddow et al., 2022); this means that
most languages in the world cannot benefit from mainstream advances and models
(Joshi et al., 2020). There is therefore a clear appeal to developing simple
architectural mechanisms for these models that are targeted at yielding
improvements in data-scarce scenarios, while interfering minimally with
mainstream preprocessing, tokenization, and training pipelines.

In the context of a low-resource language (LRL), we are often interested in
translation to and from a closely related HRL, which possible has linguistic
genealogical, regional, and cultural ties with the LRL,111This is the case,
for example, for several dialects of Arabic, all closely related to relatively
high-resource Modern Standard Arabic, and languages of the Turkic and Indic
language continua. in order to make the abundant content in HRLs available in
related LRLs. We expect that closely related languages share considerable
overlap at the subword level from cognates, borrowings and shared vocabulary
(see examples in Figure 1). Given the absence of large parallel corpora for
our language pair, we aim to leverage this shared knowledge across source and
target, intuitively, to provide “easier” routes for our MT model from source
to target sentence.

Figure 1: Translation equivalents for Bhojpuri (top) and Hindi (bottom),
demonstrating subword overlap.

Pointer Generator Networks (PGNs; See et al. (2017)) are a mechanism which
allow the model, for every output token produced, to either copy some token
from the input (“point”) or “generate” a token as per usual from the
vocabulary. PGNs have been used for a variety of problems, described in
Section 2; however, as we far as we know, this is the first work to study its
applicability to LR NMT.

We introduce a PGN mechanism into a Transformer-based NMT architecture, and
test its performance for $6$ language pairs over $4$ low-resource ranges.
Hindi-Bhojpuri (hi-bh), Spanish-Catalan (es-ca), French-Occitan (fr-oc)
represent closely-related pairs, and Hindi-Marathi (hi-mr), Spanish-English
(es-en), and French-German (fr-de) represent more distant pairs. We expect
that PGN will help most for (1) lower-resource scenarios (2) more closely-
related language pairs (3) sentence pairs with higher subword overlap. While
PGN shows improvements in certain settings, our analysis shows that they do
not come from intended sources; in particular, we find clear lack of evidence
for the above hypotheses. We discuss various factors that contribute to this
failure, highlighting fundamental challenges for LR NMT, such as noisy
datasets, mainstream tokenization practices best suited for high-resource
scenarios, as well as linguistic and orthographic complexities that may
obfuscate underlying source-target similarities. 222Code will be made public
in non-anonymous version.

##  2 Related Work

Pointer Networks were introduced to solve problems that involved permuting the
input, such as the Traveling Salesman Problem and the complex hull problem
(Vinyals et al., 2015). Their use in NLP has been largely been for monolingual
summarization, where the target may naturally contain identical spans from the
source. Cheng and Lapata (2016) present a complex hierarchical LSTM-based
model for summarization, which directly extracts sentences from the text and
words from sentences. Gulcehre et al. (2016) use pointer networks in RNN-based
sequence-to-sequence models for summarization and machine translation,
training their model explicitly to use the pointing mechnanism for uncommon
words. Gu et al. (2016) and See et al. (2017) also incorporate variants of
pointer-generator networks into RNN-based sequence-to-sequence learning for
summarization. Prabhu and Kann (2020) applied PGNs to the task of grapheme-to-
phoneme conversion via an explicit source-target mapping. Zhang et al. (2021)
proposed a pointer-disambiguator-copier (PDC) system for dictionary-enhanced
high-resource NMT, using source word translations as potential candidates for
the copying mechanism, with a disambiguator component to select appropriate
senses.

Our work is the first to examine the applicability of PGNs as facilitators in
the low-resource MT scenario, looking to exploit linguistic relationships
between the source and target in the absence of external resources. We work
with Transformer-based NMT, and make no changes to standard BPE tokenization
schemes or training objectives (unlike Gulcehre et al. (2016) and Zhang et al.
(2021)). This is so that our findings are most relevant in today’s paradigm of
generalized strategies for end-to-end multilingual MT; our mechanism can be
easily plugged into and trained with any modern (multilingual) MT pipeline.

##  3 Model

The PGN model provides two routes to the model for predicting any target
token: copying from the source or generating from the vocabulary (Prabhu and
Kann, 2020). Copy and generate distributions at step $t$ are mixed using a
learned parameter $p_{copy}^{t}$, to obtain the final probability distribution
$P^{t}$ for the target token.

$\displaystyle
p_{copy}^{t}=\sigma(\textbf{W}^{T}(\textbf{c}^{t}\oplus\textbf{d}^{t}\oplus\textbf{s}^{t})+\textbf{B})$
$\displaystyle\textbf{P}^{t}=p_{copy}^{t}\cdot\textbf{P}_{\textbf{c}}^{t}+(1-p_{copy}^{t})\cdot\textbf{P}_{\textbf{g}}^{t}$

Here, $\textbf{c}^{t}$ is the context vector, calculated as
$\textbf{c}^{t}=(\textbf{a}^{t})^{T}\textbf{e}^{t}$, where $\textbf{a}^{t}$
represents cross-attention vector, and $\textbf{e}^{t}$ contains the encoder
hidden states. $\textbf{d}^{t}$ and $\textbf{x}^{t}$ contain the decoder’s
final hidden states and input respectively, $\oplus$ denotes concatenation,
and W and B are learned weights and a bias vector respectively.
$\textbf{P}_{\textbf{c}}^{t}$ and $\textbf{P}_{\textbf{g}}^{t}$ represent the
copy and generate distributions (softmaxed logits) respectively at step $t$.
We use cross-attentions weights over source tokens for the copy logits, and
standard decoder outputs for generate logits.

##  4 Experiments

### Datasets and languages

We used the WikiMatrix (wm) corpus (Schwenk et al., 2019) for es-en, es-ca,
fr-de and fr-oc. For es-ca, we also report results on synthetic Europarl (ep)
parallel data (Koehn, 2005).333The Europarl dataset was automatically
translated into Catalan; taken from https://github.com/Softcatala/Europarl-
catalan. For hi-mr, we used the CVIT-PIB corpus (Philip et al., 2021), and for
the low-resource pair hi-bh, we use the NLLB corpus (Schwenk et al., 2021;
Team et al., 2022; Heffernan et al., 2022). See Table 1 for dataset
heuristics. Note the higher per token overlap as expected for our closely-
related group as compare to the others. The hi-bh sentences are extremely
short, and share fewer tokens than expected: in this case, this reflects badly
parallel data.444See App A for more details on datasets.

| hi-mr | hi-bh | es-en | es-ca (ep) | es-ca (wm) | fr-de | fr-oc
---|---|---|---|---|---|---|---
Avg. common tokens per sent pair | 2.51 | 1.29 | 2.81 | 6.77 | 7.17 | 1.76 | 5.06
Avg. common tokens per target token | 0.12 | 0.16 | 0.13 | 0.26 | 0.29 | 0.10 | 0.22
Avg. source sentence length | 28.54 | 6.34 | 23.43 | 26.86 | 25.38 | 19.10 | 24.09
Avg. target sentence length | 20.82 | 7.98 | 20.77 | 26.26 | 24.64 | 16.93 | 23.43
Table 1: Statistics on commons source-target tokens in our datasets. wm:
WikiMatrix, ep: Europarl.

### Experimental Settings

For all language pairs, we performed experiments on dataset subsets of $5k$,
$15k$, $30k$, and $60k$ sentences and test sets of $5000$ sentences, trained
until convergence, with tokenizer size $16000$. We computed baseline results
(NMT) on standard encoder-decoder NMT. All PGN and NMT models use $6$ encoder
and decoder layers, $4$ attention heads, and a hidden size of $512$.

##  5 Results and Discussion

### Improvement patterns

|  |  5K |  15K |  30K |  60K |  Avg. $\Delta$
---|---|---|---|---|---|---
hi-mr |  NMT |  3.4 |  7.4 |  11.9 |  16.3 |
|  PGN |  3.4 |  7.9 |  12.6 |  18.4 |  +0.8
hi-bh |  NMT |  5.6 |  6.1 |  9.7 |  15.8 |
|  PGN |  4.0 |  8.3 |  11.4 |  12.7 |  -0.2
es-en |  NMT |  9.8 |  30.4 |  38.3 |  41.7 |
|  PGN |  9.4 |  30.0 |  38.5 |  42.2 |  0.0
es-ca(wm) |  NMT |  35.4 |  50.9 |  54.3 |  56.4 |
|  PGN |  34.7 |  51.2 |  54.1 |  57.1 |  0.0
es-ca(ep) |  NMT |  62.6 |  70.6 |  73.2 |  73.6 |
|  PGN |  62.5 |  71.6 |  74.0 |  74.2 |  +0.6
fr-de |  NMT |  3.5 |  10.8 |  19.5 |  27.2 |
|  PGN |  3.7 |  11.1 |  20.1 |  27.2 |  +0.3
fr-oc |  NMT |  24.7 |  42.2 |  45.5 |  48.7 |
|  PGN |  24.8 |  43.4 |  46.4 |  48.5 |  +0.5
Table 2: spBLEU across dataset sizes (#sents). Closely-related pairs are underlined. wm: WikiMatrix, ep: Europarl. |  | hi-mr | hi-bh | es-en | es-ca | fr-de | fr-oc |  Avg. $\Delta$
---|---|---|---|---|---|---|---|---
L | NMT | 7.6 | 10.4 | 26.6 | 44.9 | 11.9 | 35.5 |
| PGN | 8.2 | 13.5 | 24.4 | 44.6 | 10.3 | 35.3 | -0.7
H | NMT | 28.9 | 22.1 | 69.8 | 72.4 | 52.2 | 65.2 |
| PGN | 29.3 | 18.3 | 69.9 | 71.6 | 51.6 | 64.3 | -1.0
Table 3: spBLEU scores on test sets with low (L) and high (H) density of
shared source-target subwords.

See results in Table 2.555We report spBLEU since our approach attempts to
benefit performance on shared subwords. Our results are not directly
comparable to those in the literature due to differences mentioned in Section
2 and the size of training bitext ($2M$ in (Gulcehre et al., 2016), $~{}1M$ in
Zhang et al. (2021) vs. our maximum resource setting of $0.06M$).666For a
rough idea: Zhang et al. (2021) report gains in MT of $1.5-2.5$ BLEU. We see
weak improvements in a majority of settings; however, counter to intuition,
PGN does not show a clear advantage for closely-related as compared to more
distant pairs, or for lower-resource settings.

### Controlled test sets

We test the motivating hypothesis that PGN models is likely to benefit
sentence pairs with higher subword overlap. We rank sentence pairs in our test
set by percentage of shared subwords in source and target, and construct test
subsets with low and high shared-subword density from the top and bottom $500$
sentences respectively. However, in Table 3, we see that in fact that PGN
performs slightly worse than NMT on both extremes, indicating that observed
benefits over the entire test set do not come from subword overlap.

### Usage of the copy mechanism

We record values of $p_{copy}$ to track the model’s usage of the copy
mechanism. While $p_{copy}$ values are relatively high777Note that it is
difficult to comment on absolute values of $p_{copy}$. The copying
distribution is normalized over the sentence length whereas generate
distributions are normalized over the vocabulary; even low values of
$p_{copy}$ will considerably affect the mixed distribution. for copied
subwords, numerals, and proper nouns, we often see that they they are also
high for seemingly random subwords.888See App C for visualizations of this
behaviour. Examples with counter-intuitively high values of $p_{copy}$:
quiero-vull (es-ca), behad-atishay (hi-mr). We also do not see a relationship
between the $p_{copy}$ value of a target token and the entropy of the cross-
attention distribution for that token.

A reasonable intuition about PGN training generalization is that in the
absence of any information, the model will default to copying, since this is
likely to do better on average than guesses over the entire vocabulary, and
that eventually, it will learn to generate language-specific subwords,
memorising the relevant strategy for given subwords in encoder states (used to
calculate $p_{copy}$ as shown in Section 3). However, our visualizations of
cross-attention and $p_{copy}$ usage throughout training show no evidence of
this generalization strategy. It’s possible that since initial cross-attention
distributions are noisy, and most subwords are not direct copies, the model is
discouraged from copying early on; it’s also possible that the model finds it
easier or trivial to encode copied source-target equivalents via the
“generate” mechanism and does not need an explicit copier, given that it must
additionally learn which subwords should be copied. In general, it appears
that the model uses the copy mechanism to encode a task that is not easily
interpretable. We discuss potential reasons for this below.

### Tokenization

In theory, the copier would learn best if the tokenizer behaved in a
morphologically principled manner.999e.g. given khaya-khalla (ate) in Hindi
and Marathi, we ideally want kha ##ya and kha ##lla. This will allow the
common stem kha to be copied over, while the language specific inflection
subwords can be generated. However, BPE tokenization generally results in
subword splits that may not reflect shared stems in word
equivalents.101010e.g. our trained tokenizer contains both propuesto (es) and
proposat (ca) instead of sharing the subword prop.

A natural idea here may be an investigation of morphologically inspired
tokenizers (Pan et al., 2020). However, we generally see inconclusive, at best
marginal, benefits of such tokenizers over BPE in modern neural MT (Macháček
et al., 2018; Domingo et al., 2019; Mielke et al., 2021), especially those
relying on unsupervised mophological segmentation, e.g., with Morfessor
(Creutz and Lagus, 2007) in the absence of morphological analysers; these
ideas have not been incorporated into mainstream tokenization strategies.

Recent work attempts to solve this general problem by looking at maximisation
of shared subwords in multilingual tokenizers (Chung et al., 2020; Zheng et
al., 2021; Liang et al., 2023); it’s possible that such strategies will
dovetail well with PGN mechanisms if widely adopted in the future.

### Linguistic complexities

While closely-related language show high (subword) vocabulary overlap, word
equivalences may be obscured by sound change and orthographic systems; if
these changes are word-internal, then even an ideal tokenizer will see
different stems/tokens in the source and target.111111e.g. vishwas-biswas (hi-
bh, sound change), website-Webseite (en-de, orthographic system) Further, we
may see that a word that has a cognate in its sister language is translated to
a non-cognate, due to semantic drift, or differences in idiom or usage norms
in the two languages, e.g. kitaab-pustak (hi-mr), resulting in non-identical
subword equivalences. These phenomena are often unpredictable and
unsystematic; even if not, they are not trivial to model into tokenization or
architectural strategies for MT.121212See App B for experiments with minor
variants of our approach dealing with pretrained encoder/decoder
initialization, tokenizer size, choice of attention head, and identical
source-target settings.

##  6 Challenges for LR NMT

Incorporating knowledge of linguistic relationships among closely related
data-imbalanced language pairs offers a natural strategy for mitigating data
scarcity in mainstream NMT between regional languages, and it is crucial to
understand the challenges in this realm. We show that while the PGN mechanism
offers an intuitive theoretical shortcut for translation between closely
related languages, its performance in practice is limited, potentially by the
combined effect of noisy real-world datasets containing non-literal
translations, the behaviour of standard tokenisers, as well as linguistic
complexities beyond simplified ideas of shared vocabulary and cognates. These
are inevitable hurdles to any project that attempts to use structural
linguistic knowledge to benefit NMT performance.

Further, we show that despite showing benefits in certain settings over the
entire test set, the PGN mechanism does not perform as expected on target
phenomena. The generalization mechanisms of blackbox Transformer-models are
not well understood and may not be easily guided by linguistic intuition: we
underline the importance of verifying that improvements are coming from the
intended places rather than good starts or extra parameters.

Finally, our analysis hints that PGN-like shortcuts may not be worth offering
in the first place: “easy” equivalences, a natural target of linguistic
interventions, may not be the bottleneck for LR NMT. Instead, it is more
likely that the true bottlenecks are handling precisely the above challenges,
i.e. non-systemic differences, one-off phenomena, and real-world noise in low-
resource conditions.

##  7 Limitations

While we show that our particular flavour of NMT incorporating PGNs does not
provide fundamental benefits for low-resource NMT, this is naturally not to
say that an improved variant of this idea would not work better. There are
several potential ways forward arising from our discussion of the reasons for
the failure of our method in Section 5: for example, using morphological
segmentation for tokenization to increase subword overlap, or using priors for
$p_{copy}$ so that it is encouraged for shared subwords. Previous work
provides different kinds of help to the copier: for example, (Gulcehre et al.,
2016) explicitly train the copier to copy unknown words with a separate
training objective. However, as we mention in Section 1, our motivation lay in
designing a simple architectural mechanism that can be easily integrated into
mainstream (multilingual) NMT pipelines to make them more capable for low-
resource MT, without requiring much additional language-pair specific
attention to training paradigms or tools such as morphological analysers and
bilingual lexicons, which are in any case of poor quality for low-resource
languages. We restrict our negative result to the scope set up by this
motivation.

Further, our results are limited to the $6$ language pairs that we
experimented on. While we simulate identical low-resource conditions for all
our language pairs, we clearly see the difference in absolute performances on
hi-mr or hi-bh as compared to the high-resource language pairs: the data for
the latter are simply of much better quality. This demonstrates the need to
experiment and present results on truly low-resource conditions: we include
hi-bh as such a language pair. Finally, this discussion is only relevant to
translation between closely-related languages that share a script (although
this is the predominant case), allowing for lexical similarity to be reflected
by shared subwords.

## References

  * Cañete et al. (2020) José Cañete, Gabriel Chaperon, Rodrigo Fuentes, Jou-Hui Ho, Hojin Kang, and Jorge Pérez. 2020.  Spanish pre-trained bert model and evaluation data.  In _PML4DC at ICLR 2020_.
  * Cheng and Lapata (2016) Jianpeng Cheng and Mirella Lapata. 2016.  Neural summarization by extracting sentences and words.  In _Proceedings of the 54th Annual Meeting of the Association for Computational Linguistics (Volume 1: Long Papers)_ , pages 484–494, Berlin, Germany. Association for Computational Linguistics.
  * Chung et al. (2020) Hyung Won Chung, Dan Garrette, Kiat Chuan Tan, and Jason Riesa. 2020.  Improving multilingual models with language-clustered vocabularies.  In _Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing (EMNLP)_ , pages 4536–4546, Online. Association for Computational Linguistics.
  * Creutz and Lagus (2007) Mathias Creutz and Krista Lagus. 2007.  Unsupervised models for morpheme segmentation and morphology learning.  _ACM Transactions on Speech and Language Processing (TSLP)_ , 4(1):1–34.
  * Devlin et al. (2018) Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova. 2018.  Bert: Pre-training of deep bidirectional transformers for language understanding.  _arXiv preprint arXiv:1810.04805_.
  * Domingo et al. (2019) Miguel Domingo, Mercedes Garcıa-Martınez, Alexandre Helle, Francisco Casacuberta, and Manuel Herranz. 2019.  How much does tokenization affect neural machine translation?
  * Forcada et al. (2011) Mikel L Forcada, Mireia Ginestí-Rosell, Jacob Nordfalk, Jim O’Regan, Sergio Ortiz-Rojas, Juan Antonio Pérez-Ortiz, Felipe Sánchez-Martínez, Gema Ramírez-Sánchez, and Francis M Tyers. 2011.  Apertium: a free/open-source platform for rule-based machine translation.  _Machine Translation_ , 25(2):127–144.
  * Gu et al. (2016) Jiatao Gu, Zhengdong Lu, Hang Li, and Victor OK Li. 2016.  Incorporating copying mechanism in sequence-to-sequence learning.  _arXiv preprint arXiv:1603.06393_.
  * Gulcehre et al. (2016) Caglar Gulcehre, Sungjin Ahn, Ramesh Nallapati, Bowen Zhou, and Yoshua Bengio. 2016.  Pointing the unknown words.  _arXiv preprint arXiv:1603.08148_.
  * Haddow et al. (2022) Barry Haddow, Rachel Bawden, Antonio Valerio Miceli Barone, Jindřich Helcl, and Alexandra Birch. 2022.  Survey of low-resource machine translation.  _Computational Linguistics_ , 48(3):673–732.
  * Heffernan et al. (2022) Kevin Heffernan, Onur Çelebi, and Holger Schwenk. 2022.  Bitext mining using distilled sentence representations for low-resource languages.
  * Joshi et al. (2020) Pratik Joshi, Sebastin Santy, Amar Budhiraja, Kalika Bali, and Monojit Choudhury. 2020.  The state and fate of linguistic diversity and inclusion in the NLP world.  In _Proceedings of the 58th Annual Meeting of the Association for Computational Linguistics_ , pages 6282–6293, Online. Association for Computational Linguistics.
  * Joshi (2022) Raviraj Joshi. 2022.  L3cube-hindbert and devbert: Pre-trained bert transformer models for devanagari based hindi and marathi languages.  _arXiv preprint arXiv:2211.11418_.
  * Khanna et al. (2021) Tanmai Khanna, Jonathan N Washington, Francis M Tyers, Sevilay Bayatlı, Daniel G Swanson, Tommi A Pirinen, Irene Tang, and Hèctor Alòs i Font. 2021.  Recent advances in apertium, a free/open-source rule-based machine translation platform for low-resource languages.  _Machine Translation_ , 35(4):475–502.
  * Kim et al. (2020) Yunsu Kim, Miguel Graça, and Hermann Ney. 2020.  When and why is unsupervised neural machine translation useless? In _Proceedings of the 22nd Annual Conference of the European Association for Machine Translation_ , pages 35–44, Lisboa, Portugal. European Association for Machine Translation.
  * Koehn (2005) Philipp Koehn. 2005.  Europarl: A parallel corpus for statistical machine translation.  In _Proceedings of Machine Translation Summit X: Papers_ , pages 79–86, Phuket, Thailand.
  * Liang et al. (2023) Davis Liang, Hila Gonen, Yuning Mao, Rui Hou, Naman Goyal, Marjan Ghazvininejad, Luke Zettlemoyer, and Madian Khabsa. 2023.  Xlm-v: Overcoming the vocabulary bottleneck in multilingual masked language models.  In _Proceedings of the 2023 Conference on Empirical Methods in Natural Language Processing_ , page 13142–13152, Singapore. Association for Computational Linguistics.
  * Macháček et al. (2018) Dominik Macháček, Jonáš Vidra, and Ondřej Bojar. 2018.  Morphological and language-agnostic word segmentation for nmt.
  * Mielke et al. (2021) Sabrina J. Mielke, Zaid Alyafeai, Elizabeth Salesky, Colin Raffel, Manan Dey, Matthias Gallé, Arun Raja, Chenglei Si, Wilson Y. Lee, Benoît Sagot, and Samson Tan. 2021.  Between words and characters: A brief history of open-vocabulary modeling and tokenization in nlp.
  * Pan et al. (2020) Yirong Pan, Xiao Li, Yating Yang, and Rui Dong. 2020.  Morphological word segmentation on agglutinative languages for neural machine translation.
  * Philip et al. (2021) Jerin Philip, Shashank Siripragada, Vinay P. Namboodiri, and C. V. Jawahar. 2021.  Revisiting low resource status of indian languages in machine translation.  In _Proceedings of the 3rd ACM India Joint International Conference on Data Science & Management of Data (8th ACM IKDD CODS & 26th COMAD)_, page 178–187, Bangalore India. ACM.
  * Prabhu and Kann (2020) Nikhil Prabhu and Katharina Kann. 2020.  Making a point: Pointer-generator transformers for disjoint vocabularies.  In _Proceedings of the 1st Conference of the Asia-Pacific Chapter of the Association for Computational Linguistics and the 10th International Joint Conference on Natural Language Processing: Student Research Workshop_ , pages 85–92.
  * Schwenk et al. (2019) Holger Schwenk, Vishrav Chaudhary, Shuo Sun, Hongyu Gong, and Francisco Guzmán. 2019.  Wikimatrix: Mining 135m parallel sentences in 1620 language pairs from wikipedia.  _arXiv preprint arXiv:1907.05791_.
  * Schwenk et al. (2021) Holger Schwenk, Guillaume Wenzek, Sergey Edunov, Edouard Grave, Armand Joulin, and Angela Fan. 2021.  CCMatrix: Mining billions of high-quality parallel sentences on the web.  In _Proceedings of the 59th Annual Meeting of the Association for Computational Linguistics and the 11th International Joint Conference on Natural Language Processing (Volume 1: Long Papers)_ , pages 6490–6500, Online. Association for Computational Linguistics.
  * See et al. (2017) Abigail See, Peter J Liu, and Christopher D Manning. 2017.  Get to the point: Summarization with pointer-generator networks.  _arXiv preprint arXiv:1704.04368_.
  * Team et al. (2022) NLLB Team, Marta R. Costa-jussà, James Cross, Onur Çelebi, Maha Elbayad, Kenneth Heafield, Kevin Heffernan, Elahe Kalbassi, Janice Lam, Daniel Licht, Jean Maillard, Anna Sun, Skyler Wang, Guillaume Wenzek, Al Youngblood, Bapi Akula, Loic Barrault, Gabriel Mejia Gonzalez, Prangthip Hansanti, John Hoffman, Semarley Jarrett, Kaushik Ram Sadagopan, Dirk Rowe, Shannon Spruit, Chau Tran, Pierre Andrews, Necip Fazil Ayan, Shruti Bhosale, Sergey Edunov, Angela Fan, Cynthia Gao, Vedanuj Goswami, Francisco Guzmán, Philipp Koehn, Alexandre Mourachko, Christophe Ropers, Safiyyah Saleem, Holger Schwenk, and Jeff Wang. 2022.  No language left behind: Scaling human-centered machine translation.
  * Tiedemann (2012) Jörg Tiedemann. 2012.  Parallel data, tools and interfaces in OPUS.  In _Proceedings of the Eighth International Conference on Language Resources and Evaluation (LREC’12)_ , pages 2214–2218, Istanbul, Turkey. European Language Resources Association (ELRA).
  * Vaswani et al. (2017) Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin. 2017.  Attention is all you need.  _Advances in neural information processing systems_ , 30.
  * Vinyals et al. (2015) Oriol Vinyals, Meire Fortunato, and Navdeep Jaitly. 2015.  Pointer networks.  _Advances in neural information processing systems_ , 28.
  * Zhang et al. (2021) Tong Zhang, Long Zhang, Wei Ye, Bo Li, Jinan Sun, Xiaoyu Zhu, Wen Zhao, and Shikun Zhang. 2021.  Point, disambiguate and copy: Incorporating bilingual dictionaries for neural machine translation.  In _Proceedings of the 59th Annual Meeting of the Association for Computational Linguistics and the 11th International Joint Conference on Natural Language Processing (Volume 1: Long Papers)_ , pages 3970–3979.
  * Zheng et al. (2021) Bo Zheng, Li Dong, Shaohan Huang, Saksham Singhal, Wanxiang Che, Ting Liu, Xia Song, and Furu Wei. 2021.  Allocating large vocabulary capacity for cross-lingual language model pre-training.  In _Proceedings of the 2021 Conference on Empirical Methods in Natural Language Processing_ , pages 3203–3215, Online and Punta Cana, Dominican Republic. Association for Computational Linguistics.

##  Appendix A Notes on Datasets

The Hindi-Marathi WikiMatrix dataset (Schwenk et al., 2021) has only 11k
sentences so we use CVIT-PIB (Philip et al., 2021) instead. The CVIT-PIB
corpus is automatically aligned using an iterative process that depends on
neural machine translation into a pivot language and filtering heuristics.
Eyeballing the data, we observe a considerable number of non-parallel or even
entirely unrelated sentences simply containing some words in common - in
general, this corpus is much more likely to contain rough paraphrases as
opposed to literal translations.

For Hindi-Bhojpuri, NLLB seems to be the only available parallel text data for
now (Tiedemann, 2012); this corpus has also been automatically crawled (Team
et al., 2022). The Hindi-Bhojpuri NLLB dataset contains extremely short
sentences as shown in Table 1; similarly to above, we observe a high level of
noise and non-parallel data. Such datasets naturally do not provide the most
favourable training conditions for the PGN models, which rely on literal
translations containing shared subwords to teach the copier; however, they are
realistic real-word conditions for truly low-resource languages, as we discuss
in Section 6.

The WikiMatrix dataset (Schwenk et al., 2019), which we use for Spanish-
English, Spanish-Catalan, French-German, and French-Occitan is automatically
aligned from Wikipedia content in these languages.

The Spanish-Catalan synthetic Europarl bitext is created by automatically
translating the Europarl dataset (Koehn, 2005) into Catalan using Apertium
(Forcada et al., 2011; Khanna et al., 2021). While this data probably contains
some noise due to MT errors and translationese, it’s the most likely of all
our datasets to contain literal, linear translations, and we include it as a
testbed for this purpose. It is a generally easier dataset - this is clearly
visible from spBLEU scores that our models achieve on it in Table 2.

##  Appendix B Minor Variations

### Pretrained encoder and decoder

We tried using a pretrained encoder and decoder at initialization of our model
and tested this for hi-mr and es-ca. For the former, the encoder and decoder
were initialized with Hindi BERT (Joshi, 2022), and for the latter, we used
Spanish BERT (Cañete et al., 2020). These pre-trained models are language-
specific instances of BERT (Devlin et al., 2018). Note that this means that we
also used the pretrained tokenizers of these models, of sizes $52000$ and
$31002$ for Hindi and Spanish respectively, that are only trained on the high-
resource source languages; this leads to very poor tokenization in the target
language. In general, this set of models take longer to converge due to their
size, and show only minor differences in performance. Another related idea is
to finetune NLLB or another multilingual MT model with an incorporated PGN; we
did not try this given the lack of encouraging results from these experiments.

### Single attention head

We also tried using only a single attention head to calculate $p_{copy}$ for
target tokens, with the motivation that it was maybe better to nudge a single
head to encode information about whether target token need to be copied, and
leaving other heads to generate, as opposed to asking all heads to do both
(which is the case when we average over heads). However, these models give
almost identical results as in Table 2.

### Smaller tokenizer

We hypothesize that using a smaller tokenizer size will force more splits per
token, increasing the chance that common stems will be reflected in shared
subwords. Accordingly, we tried a tokenizer size of $8000$ for hi-mr and es-ca
for the $15k$ and $60k$ settings; however, performance degrades slightly
(about $-1$ spBLEU on average) for both NMT and PGN approaches without
affecting the relative trend.

This is not altogether surprising: reducing tokenizer size only increases the
degree of splitting in words of a certain (lower) frequency range, rather than
affecting the number of splits for all words uniformly. More importantly,
while these hyperparameters are important to tune, statistical frequency-based
tokenizers behave inherently differently from morphologically-inspired
tokenizers, as discussed in Section 5, and it is not easy or perhaps possible
to achieve a good approximation of the latter by playing with the
hyperparameters of the former.

### Identical source and target

Finally, we also trained a Hindi-Hindi model, to remove the effects of noisy
translations and non-ideal tokenization of source and target token sequences
as discussed in Sections 4 and 5. In this setup, with $100\%$ overlap, the
models achieve high test scores ($74$ spBLEU) and converge to near-zero usage
of the copy mechanism. Clearly, the model still prefers to encode the identity
relationship using a generate mechanism.

Note that this setup is fundamentally different from our other scenarios -
when all tokens are copied, the model no longer needs the distinction between
two distinct processes (generating and copying), and therefore does not really
need to learn how to make this decision. However, it is still illustrative in
demonstrating that model generalization mechanisms, even for highly simplified
or trivial tasks, are often not intuitive or human-interpretable.

##  Appendix C Visualizations

See Figures 2, 3, 4, and 5 for visualizations of the PGN model’s cross-
attention distributions and values of $p_{copy}$ per target token on randomly
chosen source-target pairs using an early and late model training checkpoint.
We observe that the model does use the copying mechanism as intended in many
places, for common subwords (udaar in hi-mr, un in es-ca) as well as named
entities (Cour in fr-oc), common borrowings (computer in hi-bh), numbers (1970
in fr-oc) and punctuation. However, $p_{copy}$ values are also relatively high
sometimes for other seemingly random target tokens, e.g. canton-costat in fr-
oc.

Note that the hi-bh source-target sentence pairs are not in fact translations
of each other and exemplify the noise we discuss in Section 4 and App A. The
cross-attention distributions for the es-ca and fr-oc models are in general
much better defined and able to attend to appropriate tokens (in these example
visualizations as well as others that we looked at); this is a consequence of
the better quality of the data and models in these languages.

(a) Epoch 10

(b) Epoch 30

Figure 2: Model’s cross-attention distributions and $p_{copy}$ values for two
sentence pairs for es-ca(ep), $60k$ sentences

(a) Epoch 10

(b) Epoch 30

Figure 3: Model’s cross-attention distributions and $p_{copy}$ values for two
sentence pairs for fr-oc, $60k$ sentences

(a) Epoch 10

(b) Epoch 25

Figure 4: Model’s cross-attention distributions and $p_{copy}$ values for two
sentence pairs for hi-mr, $60k$ sentences

(a) Epoch 10

(b) Epoch 30

Figure 5: Model’s cross-attention distributions and $p_{copy}$ values for two
sentence pairs for hi-bh, $60k$ sentences


# Using Analytics on Student Created Data to Content Validate Pedagogical
Tools

###### Abstract

Conceptual and simulation models can function as useful pedagogical tools,
however it is important to categorize different outcomes when evaluating them
in order to more meaningfully interpret results. VERA is a ecology-based
conceptual modeling software that enables users to simulate interactions
between biotics and abiotics in an ecosystem, allowing users to form and then
verify hypothesis through observing a time series of the species populations.
In this paper, we classify this time series into common patterns found in the
domain of ecological modeling through two methods, hierarchical clustering and
curve fitting, illustrating a general methodology for showing content validity
when combining different pedagogical tools. When applied to a diverse sample
of 263 models containing 971 time series collected from three different VERA
user categories: a Georgia Tech (GATECH), North Georgia Technical College
(NGTC), and “Self Directed Learners”, results showed agreement between both
classification methods on 89.38% of the sample curves in the test set. This
serves as a good indication that our methodology for determining content
validity was successful.

John Kos<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta Georgia
Kenneth Eaton<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta
Georgia
Sareen Zhang<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta
Georgia
Rahul Dass<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta
Georgia
Stephen Buckley<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta
Georgia
Sungeun An<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta
Georgia
Ashok Goel<EMAIL_ADDRESS>Georgia Institute of Technology, Atlanta
Georgia

##  1 Introduction

Students in the studies of ecology and biology often need tools to make sense
of the complex system found in real-world ecology. These pedagogical tools,
often necessary for learning, simplify real-world phenomena and make it easier
to learn and understand the myriad of factors affecting real-world
circumstances can be difficult to capture with empirical approaches. Thus,
modeling ecological systems may be a more tractable approach [ENB12].

Among the tools used in ecological education are mathematical models,
conceptual models, and simulations. All three of these tools offer different
insight to help students understand complex ecological phenomena. Mathematical
models are especially useful in describing individual population dynamics that
delineate the growth or decline of species’ population over time, given
various factors such as birth rate, death rate, etc. This description of
species’ population over time often gets categorized into common curves in
order to relate them to phenomena often found in nature [LLE08]. However, when
engaging in systems thinking, mathematical models become exponentially more
complex and challenging to understand [Pie81].

Conceptual models represent the behavior of a system by defining its
components, relations, and processes [ABH+18]. This allows students, through
predictive inquiry, to formulate ”what if” hypotheses that can be formalized
with a conceptual model. Although conceptual models provide an overview of a
system and formulate a hypothesis, they cannot evaluate how those dynamics
affect the system over time. Simulation models provide a solution to this by
using complex algorithms to model ecosystem dynamics and create the most
accurate predictions [20222]. In short, simulations provide the ability to
evaluate predictive inquiry.

The Virtual Ecological Research Assistant (VERA) [oTDL23] is an ecological
modeling software that allows for an end-to-end evaluation process from
defining the conceptual model, see Fig. 1, to observing the simulation
outputs, see Fig. 2 [ABH+18]. This provides a near-instantaneous loop of
hypothesis formation and evaluation. As outputs, VERA displays changes in the
species’ populations over time. In Section 2.1, we describe how the VERA
application works in more detail.

Figure 1: An example of a predator-prey conceptual model in VERA. Figure 2: An
example of the simulation output from a conceptual model in VERA

Simulations for pedagogical purposes are an increasingly popular method for
teaching students across disciplines, including but not limited to social-
ecology, epidemiology, and economics [JMA14, DB15, DlTOC+21]. Despite its
popularity, there is a gap in educational simulation systems research with
concern to assessing the internal content validity of the simulation based
educational tools [CH16]. By corroborating the expected behavior of
encountering a known group of curves found in ecological modeling with the
unsupervised bottom up approach of hierarchical clustering, we are proposing a
methodology for validating the use of simulations used in education, with an
emphasis primarily on domains with a heavy focus on modeling.

This methodology is also particularly useful for agent-based simulation
systems for which the complexity of the model may cause outcomes unaligned
with the educational purpose. This is due to the fact that the agent based
system in VERA is attempting to simulate the complex system of ecology.
Complex systems are defined by Ottino [Ott03] as “a system with a large number
of elements, building blocks or agents, capable of exchanging stimuli with one
another and with their environment”. Interactions within complex systems
create emergent patterns. Although VERA takes many fundamental principles from
ecology, that does not necessarily guarantee that the emergent patterns will
be the same. Ottino [Ott03] continues “The common characteristic of all
complex systems is that they display organization without any external
organizing principle being applied.” The convergence of expected outcomes
between modeling–based disciplines and agent–based simulations of those
disciplines can serve as a validation of the models used, especially in cases
where analysis of the mathematical models reaches a point of infeasible
difficulty [HCR06].

In this paper, we demonstrate a methodology by which the content validity of
pedagogical tools, specifically interactive modeling environments in domains
that deal heavily with modeling, curves and simulation, can be justified. As
an example, we apply this methodology to VERA due to its function as a tool
that can produce simulation output time series that match common population
curves that students learn about in biology ecology classes. By using VERA,
students often create conceptual models and simulate them in order to get a
time series that aligns with what would be the population graph of a
mathematical model. In short, we will demonstrate that the conceptual model
and simulation created by VERA will also reaffirm students’ understanding of
mathematical ecological models. To do this, we implement two approaches,
Hierarchical Clustering and Curve Fitting. Hierarchical clustering functions
as a bottom-up unsupervised method meaning that natural patterns in the data
will emerge[MC17, Kas17]. Curve Fitting was chosen because it is a supervised,
top-down approach to classifying the time series. This method means that we
can take inspiration from known curves used in ecological modeling. By finding
a series of patterns in the hierarchical clustering that match, to a necessary
threshold of accuracy, a group of chosen curves used in ecological modeling we
can affirm that the simulation output of VERA would be useful for reinforcing
population curves that a student learns in the classroom.

The VERA output of population graphs in the form of time series data can be
difficult to interpret due to the presence of noise, feature correlation, and
high dimensionality due to the length of the time series [KTGA13, AMST11]. As
mentioned earlier, time series analysis introduces unique challenges compared
to other forms of data. One of the more significant factors is that VERA
models and simulates complex ecological systems, so therefore outcomes may be
vastly different and there may be potential noise among relatively consistent
outcomes [AGG+10]. Another important consideration is that, in general, data
samples may consist of different lengths of time, and different levels of
starting populations for the species. When combined, these and other factors
can result in great difficulty adapting to the variability of samples that
should belong to the same class [BA11, LVVC11].

To address these issues, we develop and compare the performance of two
classification methodologies: (1) hierarchical clustering with dynamic time
warping and (2) curve fitting, to identify patterns in species population
outputs. We first explore the use of hierarchical clustering with dynamic time
warping to extract features in order to identify patterns in the population
outputs [Łuc16]. Hierarchical clustering seeks to identify groups of similar
data samples in unlabeled data using a distance/similarity measure [MC17].
Hierarchical clustering can be agglomerative, meaning it starts with every
sample in their own cluster and merges them until they form one cluster or
divisive, which does the opposite [OES07]. Łuczak explored using hierarchical
approaches on a variety of time series data sets [Łuc16].

Another popular approach that can be used to classify time series data is
curve fitting, which consists of applying mathematical functions to fit the
data to a particular curve. These can be easier to adapt than empirical
methods and can effectively suppress noise [ZWX+20]. However, they can still
struggle to adapt to data that does not align with standard function curves
[ZWX+20]. Curve fitting has been used extensively for extracting vegetation
phenological metrics with different curve types including Logistic [CCST15]
and Gaussian [JE02] curves. The curve-fitting approach took inspiration from
common curves used to describe population dynamics found in ecological
modeling. As part of our initial experiments, see Section 3, we considered six
patterns, exponential growth, capped growth, oscillation, dying off, Gaussian
and constant. A wide variety of possible curves were tested, however, the set
of curves with the greatest success rate were those that match the patterns
described above.

Across our tested sample, the curve fitting labeling and the hierarchical
clustering with a second KNN corroboration step converged on the same curve
type 89.38 percent of the time in the test set. For the rest of this paper, we
will refer to this convergence between methods as accuracy.

##  2 Background

###  2.1 VERA

VERA is an ecological modeling framework that allows users to create
conceptual models and simulate their outcome [ABH+18]. The conceptual models
are expressed in the Component-Mechanism-Phenomenon (CMP) language, which is
an extension of Structure-Behavior-Function models. [JGP14, GJ15, GRV09]. This
format of modeling allows for components such as biotics, abiotics, and
habitats to be defined as well as relationships between them [ABH+18]. Each
component type has a set of adjustable parameters allowing for specialization
to particular species. In order to assist in using valid parameters, users are
able to lookup species with preset values from the Encyclopedia of Life
(EOL)[PWL+14]. Once a model has been built, VERA translates the conceptual
model into Netlogo [TW04b] primitives using a compiler. Running the simulation
allows the user to see a plot see Fig. 2 of the population changes over time
and they can then iterate through model-simulate-refine loops [WF90] to refine
their model.

VERA is used by a variety of stakeholders including students from Biology and
Computer Science students at Georgia Institute of Technology (GATECH), Ecology
Students at North Georgia Technical College (NGTC), and self directed learners
that can find the site through the Smithsonian run EOL. The variety in
stakeholders corresponds to learners with different goals. For example, the
Computer Science students at Georgia Tech may use VERA to explore a scientific
method of thinking as part of a cognitive science class, while an ecology
student at North Georgia Technical College may model in VERA to explore the
possibility space of an ecology. Self-directed learners are different from
both aforementioned groups in the facts that their goals remain unknown.
Nevertheless, the model-simulate-refine loops implicit in VERA’s design
encourage learning where feedback is necessary for understanding the
conceptual modeling system, and therefore the effects of ones prospective
modeled hypothesis.

###  2.2 Related Works

To the best of the authors’ knowledge, we have not come across studies that
categorize the time series of ecological simulation models in the context of
education and its value in the use of pedagogical tools. This research lies at
an understudied convergence of different topics. As described by Cook [CH16],
while simulation based tools are becoming increasingly more popular, there is
a lack of validation work proving these tools usefulness. Since the publishing
of the above source, there has been a number of studies that discuss the
validity of simulation tools in educational contexts, but primarily in the
medical field, and many of those leave validation studies for future work
[BSM+21, MTD+18, KW21]. However this does not mean that there is not value in
doing validity studies for educational tools in other fields.

Borrowing from a validity framework developed by Messick used previously for
psychological construct assessment, but cited in the studies above, there are
six aspects of validity: content, substantive, structural, generalizability,
external, and consequential [Mes95]. This research particularly focuses on
content validity which is concerned with whether the assessment, or in this
case the interactive modeling environment, resembles the domain in questions.
As shown by Cook in an earlier paper concerning what counts as validity
evidence, the two most popular forms of content evidence were ”group consensus
or expert review”, and creating an instrument ”based on (or modified from) a
previously validated instrument” [CZH+14]. This research only concerns the
second, however the first is mentioned briefly in discussion. The previously
validated instrument in this case is the domain of ecological modeling, and by
reproducing validity tested methods in VERA we hope to extend content validity
from the methods to the software itself.

Approaching the problem from the opposite direction, Holcombe suggests that
”formally defined agents [within agent based systems] could be a basis for the
validation and verification of models in the longer term” [HCR06]. Framed
around the VERA software, this would take the form of ecological modeling
validating work in their field using a tool such as VERA. There are a number
of studies in environmental studies and ecology that take directly this
approach with VERA’s underlying application NetLogo [JRKB19, TKG12]. Other
existing research with this focus exists in the fields of economics, and
social sciences [WFM07, OR06]. However, a literature review did not reveal any
research with a focus on validation for pedagogical purposes in non-medical
domains. Rather than attempting to validate known models through agent based
simulations, this research seeks to use the convergence of both analytics
tools as a method to determine pedagogical content validity in interactive
learning environments.

###  2.3 Ecological Population Behaviors

Curves are used to describe population behaviors in the field of ecological
modeling. The most foundational of these curves is the exponential growth
curve which describes a closed population, modeled continuously with a
positive instantaneous rate of increase [Got95a]. Exponential growth can
typically be thought of as a new species being introduced to an environment,
like an invasive species. This species may have ample food supplies and no
natural predators initially, leading to rapid population growth [SB13].
Certain types of living organisms can also be geared more towards exponential
growth such as bacteria because of their quick reproductive cycle.

Conversely, a negative instantaneous rate of increase creates the exponential
decay curve. This occurs when a species is exposed to an environment that does
not support it, the population can rapidly decline [WB97].

However, the exponential curves assume infinite resources for consumption,
unlimited reproduction, and density independent birth and death rates.
Limiting those assumptions by capping resources and modifying the birth and
death rates to be density dependent creates the foundation for the logistic
curve [Got95b]. Another curve which is used to describe the slowing growth
rate of populations over time is the sigmoid curve. In ecological modeling,
there is discussion of whether the the sigmoid curve would be a better
representation of the circumstances used to describe the logistic curve above
[Gam98]. The sigmoid curve is also described by an exponential curve which
slows in growth as it reaches a population cap. For the purposes of this
paper, when describing curves, we will refer to both the logarithmic curve and
the sigmoid curve as capped growth as that is an accurate descriptor for the
circumstance which creates both curves.

Another model that includes species density and non closed populations in the
form of environmental gradients is the Gaussian curve [JHC74]. A Gaussian
curve could represent a species moving into a new environment where resources
are plentiful and then consuming at an unsustainable rate such that the
population then falls to a more sustainable level.

Populations in nature are also subject to temporality in that there may be
resource abundant seasons or migration and location dependent mating. In the
case of resource abundant seasons, the consumer can be said to be CR coupled
with their resource [J.04]. This creates our oscillation curve. In the case of
VERA, our oscillation curve is caused by the interactions with predator
populations which becomes cyclical when the random interaction between
populations reaches a stable equilibrium such that neither population falls
below a critical minimum population size. These behaviors explain five of the
six patterns we identify as having relevant ecological meaning as these five
curves deal primarily with the biotic component of VERA. We are also proposing
another curve, constant, for use in classification as it accounts for the
abiotic component in VERA. [TW04a]. This could be average number of square
inches of rain per year or vital non-biotic nutrients in the soil. Lastly we
will be using a category for all population graphs which do not fit into the
above curves for which we will refer to as outliers.

##  3 Methodology

###  3.1 Data processing

The data set was manually downloaded from VERA, first by identifying a user,
cloning each of their models and then running a simulation to create the data.
The time series created by the simulation was then downloaded as a CSV.

In order to account for variations in the population scales in the graphs,
which can go from 0-25,000, we scale the data from zero to one based on the
max population value in the graph. All data lengths are also cut to 400 values
(400 months in simulation), so that the samples are all the same length. This
is the default simulation length for VERA, however students are able to
lengthen or shorten the default simulation length. Graphs shorter than 400
data points were not included in this sample.

Three different user sets were looked at in this analysis, and they were
chosen primarily based on the reasoning that their differing motivation for
using VERA would lead to differing types of modeling behavior. The first user
set was Computer Science students from Georgia Tech (GATECH). As part of the
cognitive science course, students are asked to use VERA and describe how it
can assist with ”scientific thinking”. Twenty students were then randomly
selected from the course rosters of the Spring and Summer 2022 sessions for
the course. Every single model from a chosen student was then used in the data
set resulting in 177 models and 685 species/curves total.

The second user set was students studying ecology at North Georgia Technical
College (NGTC). Every student who took an ecology course that used VERA across
the Summer and Fall 2022 was included in the data set. In total this was 19
students, creating 57 models, including 194 species.

This last user group was Self-Directed Learners (SDLs). This is our
classification for learners who are not tied to any known classroom use of
VERA. In order to determine this, all users of VERA were pulled, and any that
contained a GATECH email, a NGTC email, or had email that had a first and last
name that matched a known student were removed from the data. This was
evaluated by any new users from January 1st of 2022 to September 1st 2023. The
resulted in 55 users, of which only 23 had actually created a model. From
those 23 users, 29 models were pulled, containing 92 species.

A superset of all of the above user groups was also evaluated consisting of
263 models and 971 total species. Below is a table showing the curves broken
up by type and user group.

| GATECH | NGTC | SDLs | All
---|---|---|---|---
Curve | # | % | # | % | # | % | # | %
Exponential | 17 | 2.48 | 11 | 5.67 | 1 | 1.08 | 29 | 2.99
Capped Growth | 48 | 7.01 | 18 | 9.28 | 9 | 9.78 | 75 | 7.72
Exponential Dying | 287 | 41.89 | 92 | 47.42 | 41 | 44.56 | 420 | 43.33
Oscillation | 28 | 4.09 | 28 | 14.43 | 13 | 14.13 | 69 | 7.10
Constant | 145 | 21.17 | 2 | 1.03 | 15 | 16.30 | 162 | 16.68
Gaussian | 127 | 18.54 | 33 | 17.01 | 9 | 9.78 | 169 | 17.40
Outlier | 33 | 4.82 | 11 | 5.67 | 4 | 4.34 | 47 | 4.84
Total | 685 | 100 | 194 | 100 | 92 | 100 | 971 | 100
Table 1: Graphs broken up by curve type and data set, number # and percent %

###  3.2 Old Training Set

The lab had previously done research into self directed learners in order to
better understand the types of models created by people outside a classroom
setting where we know the models that they will create due to the assignment.
Overall, it was discovered that many users copied the exemplar models, or
created models exhibiting the same basic structure but used their own species
. [S.] To understand the behaviors of the models, we wanted to align the
simulation behaviors with what is expected for certain ecological
relationships. For example, oscillation may represent a predator-prey
relationship where the populations cycle with each other. While these types of
relationships can be easily defined in an ecological sense, it was not clear
how to align them with the simulation outputs. In order to do this, we
explored both top down and bottom up clustering of the simulation outputs to
search for meaningful patterns. While searching, we utilized the silhouette
score from KNN to identify meaningful clusters. After some manual work
comparing the strongest clusters to our ecological behaviors, we identified
six patterns present in the different clusters. These are the six categories
we use: constant, dying off, oscillation, exponential growth, and capped
growth, and gaussian. The curve fitting approach was then developed with these
meaningful clusters in mind, however many different other types of curves were
tested but were found to be unsuccessful. In the figure below, we can see all
the different curve types, not including constant which would just be a
straight line across the top. While this classification list is non-exhaustive
of all the curves in ecological modeling, we will see that it is comprehensive
enough for testing purposes.

Figure 3: Curve types as labeled (1) dying, (2) capped growth, (3)
exponential, (4) oscillation, (5) gaussian, (6) outlier

The old training set was comprised of self directed learners from January 1st
of 2019 to December 31st of 2021. It was made up of 197 models, which
comprised 724 species. We ran the two methods on this data set originally and
use it to set the hyperparameters for the system.

###  3.3 Hierarchical clustering

In order to cluster the data, we use agglomerative hierarchical clustering
with complete linkage, meaning clusters are merged based on the farthest
distance between points in the clusters. Similar to work by An, we use dynamic
time warping with euclidean distance in order to calculate a distance matrix
on the training data. [S.] To isolate the patterns that are subject variation,
we remove all of the constant curves from the clustering. This is done by
creating a set in python with all the values from the curve. If the set length
is one, then we know that the values of the curve do not change and therefore
it is a constant curve.

Additionally, the data set is broken into a training set made up of 70% of the
data determined by hierarchical clustering, and a test set made up of the
remaining 30% determined by KNN on the clusters. The purpose of this is two
fold. The first is to evaluate the quality of the hierarchical clusters. If we
introduce new curves and run KNN to determine their similarity to existing
clusters, then the introduced curve should fall into a similar cluster to that
of their curve type. Second, by evaluating an introduced curve by KNN of the
medoid curve, a representative curve which has the least distance between
itself and the other curves, in each cluster, we can speed up the process of
comparison, allowing for quicker classification in the future. The KNN process
is described further in the section below.

With the remaining curves, a dendrogram is created using the scipy dendrogram
function and then it is flattened using the scipy fcluster function. This
separates the dendorgram into clusters based on euclidean distance. Within
each cluster a prototypical representative, the mediod, is found determining
which curve has the least distance between itself and all of the other
members. For classification, clusters are categorized based on the label
attached to the representative curve based on the label from the curve fitting
algorithm. The label for each member of the cluster is then also evaluated in
order to see if it matched the representatives’ label in more than 55% of
cases. If the 55% threshold was not met, the entire cluster was classified as
outlier.

Figure 4: Hierarchical Clustering Process

####  3.3.1 KNN Classification

Since clustering is an unsupervised problem, we needed to develop another
approach in order to perform classification on new samples using the clusters.
For this, we first chose to pick medoid samples from each cluster, this allows
us to speed up the computation time by only comparing against a subset of the
data, since we desire to use this analysis for live feedback in VERA. Medoids
are chosen as the data point with the minimum average distance to all other
samples in a cluster. Given a new data sample, we calculate the distance from
this sample to all of the medoids and assign the class as the same as the
medoid it has a minimum distance to. In this sense, the classification can be
thought of as performing KNN on the medoids where k=1. There is also a minimum
distance threshold that identifies a sample as an outlier if it is not close
enough to any of the medoids. This minimum distance threshold value was set to
5.

###  3.4 Curve Fitting

Curve fitting was undertaken by using the expected fundamental curves and
inserting their mathematical equation, and the time series into Scipy’s curve
fit function. Curve fit then tries to return parameters for the mathematical
equation that minimizes the residual sum of squares such that it fits to the
given time series. Each time series was then tested against every curve type
in order to detect which curve type minimizes error.

Figure 5: Visualization of Curve Fitting

Due to the variability in scale of the population graphs and the
multidimensionality of the parameter space, four methods were taken to ensure
that the function was successful in finding matching curves. First, the
population graphs were normalized to the maximum value of the graph. This is
valuable because the curve fit function searches in depth for any parameter
values between integers, and this shrinks the possibility space. Limiting the
search such that the range is zero to one instead of zero to twenty five
thousand, which is around the largest maximum value in the VERA system, means
that when curve fit searches into multiple decimal places for a given
parameter, it is more likely to search successfully. Second was offering a
prototypical set of variables for each curve. This limited the bounds such
that the higher probability for success parameter space was searched first.
The prototypical set of variables, for example, meant that sin waves with a
period exponentially longer than the time series length were not searched.
Third the parameter space was divided and the function was executed for each
subsection of the parameter space. This emphasized the divided parameters
allowing for a more fine grain search as more search cycles were spent in a
smaller area. Lastly, for curves where simple rules properly outline their
behavior, a rule based system was used. In practice, this means that the
constant and exponential dying curves were evaluated using rules. Constant was
evaluated by checking if all values in the time series were the same and not
zero, and exponential dying was evaluated by checking if the population
reached zero at any point and remained there, or close, for the remainder of
the time series. Closeness to zero was determined by the hyperparameter of 4%
of the maximum population which was determined during the old training set.

In order to detect outliers, error was evaluated by measuring the residual sum
of square between the fit curve and the given population graph. A
hyperparameter of an error threshold was then manually tuned on the old
training set in order find best results, which gave a value of 5.7. Any curve
above that threshold was then classified as a outlier.

##  4 Results

Preliminary findings in the form of the silhouette scores are discussed
briefly followed by the outcomes of the this research. Our major results
includes a description of the results as well as an analysis of each of the
user sets and superset. The analysis will consist of the dendrograms for each
set as well as a description for both of the classification steps.

Figure 6: Dendrogram for all datasets. As numbered (1) superset, (2) GATECH,
(3) NGTC, (4) SDLs.

###  4.1 Silhouette Plot

Below are the silhouette scores for Kmeans clusters on the test dataset. This
helped us decide on our different curve types. As seen in the table below 7
clusters is the max number before a large drop off in score. This score then
recovers in higher N values, but meaningful ecological differences between the
clusters were not differentiable. These represent our six curves plus
outliers.

Clusters | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10
---|---|---|---|---|---|---|---|---|---
Silhouette Score | 0.7292 | 0.6319 | 0.6502 | 0.6363 | 0.6478 | 0.6471 | 0.6191 | 0.6436 | 0.6444
Table 2: Silhouette scores based on cluster number (higher better)

In the figures below we can see a representation of our clusters. For the
seven clusters as seen on the left, we really only have one highly negative
cluster, cluster 5. This represents our outliers which are not easily
classified. As we move up to 8 clusters, we develop another highly negative
cluster in cluster 7 which does not add any information value.

Figure 7: A comparison between 7 and 8 clusters

###  4.2 Accuracy Findings

User Set | GATECH | NGTC | SDLs | ALL
---|---|---|---|---
Method Accuracy (%) | 88.43 | 71.65 | 70.65 | 89.38
Table 3: Table showing user sets and the percentage of agreement between the
two methods on test set

As shown in Table 3, the accuracy decreases when clustering is applied to the
smaller sets of data. This is somewhat expected as more data would make for
more numerous clusters as well as more defined clusters.

Using 30 as the distance measure for which the clusters are evaluated, as
decided during the old training set, a total of 66 clusters were created. Of
these clusters 9 were labeled capped growth, with 59 of 61 of the curves being
labeled as capped growth. 4 of the clusters were exponential growth, with 17
of the 17 curves in the clusters being labeled as such. 11 of the clusters
were labeled exponential dying with 309 of the 320 curves matching.
Oscillation made up 11 clusters with 57 of the 71 curves matching. Gaussian
made up 22 of the clusters with 124 of the 136 curves matching. Outliers made
up 9 of the clusters with 28 of the 32 matching. This method of classification
had a success rate of 93.24%

10 additional of the above clusters were classified as outliers due to not
meeting the 55% matching threshold per cluster. These were clusters 3, 17, 21,
27, 30, 33, 40, 53, 54, 65. Of these, 7 were originally Gaussian clusters, 1
was a dying cluster and 1 was a Oscillation cluster. For the KNN validation
test set, the accuracy fell to 89.38%. This is followed by a confusion matrix,
where the rows are the expected label and the columns are the final
classification to the right of the total, and then the number of correct and
incorrect to the left of the total.

Label / Class | Out | Exp | Cap | Die | Osc | Gau | Con | Total | # Correct | # Incorrect
---|---|---|---|---|---|---|---|---|---|---
Outlier | 3 | 2 | 6 | 1 | 0 | 0 | 0 | 12 | 3 | 9
Exponential | 0 | 8 | 0 | 1 | 0 | 1 | 0 | 10 | 8 | 1
Capped Growth | 0 | 0 | 22 | 0 | 0 | 0 | 1 | 23 | 22 | 1
Dying | 0 | 0 | 0 | 117 | 0 | 3 | 0 | 120 | 117 | 3
Oscillation | 0 | 0 | 0 | 0 | 15 | 4 | 0 | 19 | 15 | 4
Gaussian | 0 | 1 | 1 | 2 | 4 | 42 | 0 | 50 | 42 | 8
Constant | 0 | 0 | 0 | 0 | 0 | 0 | 54 | 54 | 54 | 0
Table 4: Confusion matrix for test set of superset, with number of correct and
incorrect attached.

Analyzing the Georgia Tech (GATECH) data set, a total of 35 clusters were
created. Of these clusters, 5 were labeled capped growth, with 39 of 48 of the
curves being labeled as capped growth. 3 of the clusters were exponential
growth, with 13 of the 22 curves in the clusters being labeled as such. 8 of
the clusters were labeled exponential dying with 209 of the 216 curves
matching. Oscillation made up 2 clusters with 8 of the 14 curves matching.
Gaussian made up 13 of the clusters with 86 of the 97 curves matching.
Outliers made up 5 of the clusters with 20 of the 23 matching. This method of
classification had a success rate of 90.57%

7 additional of the above clusters were classified as outliers due to not
meeting the matching threshold per cluster. These were clusters 4, 6, 7, 9,
20, 25, 27. Of these, 3 were exponential, 1 was dying, 2 were Gaussian, and 1
was oscillation. For the KNN validation test set, this fell to 88.46%.

Label / Class | Out | Exp | Cap | Die | Osc | Gau | Con | Total | # Correct | # Incorrect
---|---|---|---|---|---|---|---|---|---|---
Outlier | 2 | 0 | 0 | 0 | 1 | 1 | 0 | 4 | 2 | 2
Exponential | 3 | 0 | 1 | 0 | 0 | 0 | 0 | 4 | 0 | 4
Capped Growth | 1 | 0 | 12 | 0 | 0 | 0 | 0 | 13 | 12 | 1
Dying | 0 | 0 | 0 | 84 | 0 | 1 | 0 | 85 | 84 | 1
Oscillation | 2 | 0 | 0 | 0 | 1 | 7 | 0 | 10 | 1 | 9
Gaussian | 2 | 0 | 1 | 3 | 3 | 40 | 0 | 46 | 40 | 6
Constant | 0 | 0 | 0 | 0 | 0 | 0 | 44 | 44 | 44 | 0
Table 5: Confusion Matrix for test set of GATECH with number of correct and
incorrect attached.

Analyzing the North Georgia Technical College (NGTC) data set, a total of 11
clusters were created. Of these clusters four were labeled capped growth, with
19 of 22 of the curves being labeled as capped growth. One of the clusters was
exponential growth, with 6 of the 8 curves in the clusters being labeled as
such. Two of the clusters were labeled exponential dying with 66 of the 85
curves matching. Oscillation made up one cluster with 8 of the 9 curves
matching. Gaussian made up three of the clusters with 11 of the 16 curves
matching. Outliers made up one of the clusters with 4 of the 5 matching. This
method of classification had a success rate of 78.62%

1 additional of the above clusters was classified as an outlier due to not
meeting the matching threshold. This cluster was cluster number 4, which was a
Gaussian cluster. For the KNN validation test set, this fell to 66.10%.

Label / Class | Out | Exp | Cap | Die | Osc | Gau | Con | Total | # Correct | # Incorrect
---|---|---|---|---|---|---|---|---|---|---
Outlier | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0
Exponential | 0 | 2 | 0 | 1 | 0 | 0 | 0 | 3 | 2 | 1
Capped Growth | 0 | 0 | 3 | 0 | 0 | 0 | 0 | 3 | 3 | 0
Dying | 0 | 0 | 0 | 28 | 0 | 0 | 0 | 28 | 28 | 0
Oscillation | 5 | 0 | 0 | 0 | 3 | 0 | 2 | 10 | 3 | 7
Gaussian | 3 | 0 | 0 | 1 | 1 | 1 | 7 | 15 | 1 | 14
Constant | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0
Table 6: Confusion Matrix for the test set of NGTC with number of correct and
incorrect attached.

Analyzing the Self Directed Learners data set, a total of 5 clusters were
created. Of these clusters 2 were labeled capped growth, with 8 of 8 of the
curves being labeled as capped growth. 0 of the clusters were labeled as
capped growth. 2 of the clusters were labeled exponential dying with 31 of the
41 curves matching. Oscillation made up 1 cluster with 5 of the 7 curves
matching. Gaussian made up 1 of the clusters with 2 of the 3 curves matching.
Outliers made up 0 of the clusters. This method of classification had a
success rate of 81.35%

2 additional of the above clusters were classified as an outlier due to not
meeting the matching threshold. These were clusters 1 and 3, which were dying
and Gaussian respectively. After the KNN validation test set, this fell to
64.29%.

Label / Class | Out | Exp | Cap | Die | Osc | Gau | Con | Total | # Correct | # Incorrect
---|---|---|---|---|---|---|---|---|---|---
Outlier | 0 | 0 | 1 | 0 | 0 | 0 | 0 | 1 | 0 | 1
Exponential | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0
Capped Growth | 1 | 0 | 2 | 0 | 0 | 0 | 0 | 3 | 2 | 1
Dying | 0 | 0 | 0 | 8 | 0 | 0 | 2 | 10 | 8 | 2
Oscillation | 2 | 0 | 0 | 0 | 4 | 0 | 1 | 7 | 4 | 3
Gaussian | 2 | 0 | 0 | 0 | 0 | 0 | 1 | 3 | 0 | 3
Constant | 0 | 0 | 0 | 0 | 0 | 0 | 4 | 4 | 4 | 0
Table 7: Confusion Matrix for test set of SDLs with number of correct and
incorrect attached.

##  5 Discussion and Future Work

###  5.1 Accuracy of different datasets

As seen in the results section, the accuracy of the datasets is highly
correlated with the size of the dataset. In smaller datasets, such as North
Georgia Technical College and Self Directed Learners, the accuracy loss from
the Clustering step to the KNN validation step was the largest, 12.52% and
17.06% respectively. The small cluster size meant that many clusters in these
datasets became infeasible during the KNN validation step. This is seen in the
Oscillation and Gaussian categories for the North Georgia Technical College
dataset and the Exponential and Gaussian curves in the Self Directed Learners
dataset. This even occurred in the second largest data of Georgia Tech, where
all the exponential curves were classified incorrectly.

In the largest dataset, we began to see a convergence of the clustering and
KNN accuracies. The superset of data had an accuracy of 93.24% during the
clustering training step and then an accuracy of 89.38% during the KNN
validation step. This is only a 3.86% difference, which is significantly
smaller than the smaller datasets. This indicates that if we continue to
increase our dataset size, we should see both values approach a similar
accuracy percentage.

The largest challenge going forward, and a possible avenue for future work is
the classification of outlier curves. The increased dataset size also meant
that fewer clusters were classified as outliers during the step where 55% of
curves in a cluster had to share a label with the medoid curve or else it was
classified as an outlier. 15.15% of clusters were classified as outliers
during this step in the superset, 20% percent were classified as outliers for
Georgia Tech, 9% were classified as outliers in the for North Georgia
Technical College (NGTC), and 40% for self directed learners. The lower
percentage for NGTC outlier is unsurprising due to the level of guidance that
was given while creating VERA models, their assignment is much more closed in
scope. Nevertheless, a larger dataset should converge on specific clusters for
outlier curves of similar types, while simultaneously classify less clusters
as outlier due to a lack of confluence in cluster labels. Singular cases of
outlier curves that do not follow a currently unknown outlier pattern will
still be caught during the KNN step due to its likely distance from other
curve types.

###  5.2 Education Impacts

This paper’s analysis of VERA’s simulation output is educationally valuable in
three ways. Firstly, the clustering demonstrates a mapping between VERA’s
simulation outputs and common ecological curves that some students may be
learning in the classroom, validating VERA’s use. Secondly, and by extension,
the output classes can be directly presented to learners to help them better
understand the behavior of their models. For users unfamiliar with common
ecological curves, the output classes could give them a vocabulary stepping
stone to better understand ecology. Lastly, the largest motivational purpose
of this work was in order to integrate the simulation output analysis into
VERA’s personalized coaching system. This would function by analyzing the
user’s model in comparison to other models in the database for which the
simulation output is already known. The comparison would then allow us to
suggest changes to the model which may assist the user in their learning
process. This could be useful in situations where a user is changing
parameters or relationships in the conceptual model but still seeing the same
class of simulation output, and therefore not gaining any new knowledge about
the effects of their conceptual changes. For example, when a user may have a
model that consistently has species that are dying out, the personalized
tutoring system may then suggest a node to be added or a parameter to be
changed in order for the simulation to output an oscillation. The change in
simulation output, especially if the curves are labelled and viewable by the
student, gives the student the ability to better conceptualize what changes in
the ecology caused the change in the output.

###  5.3 Validation

For students who may be using VERA in biology or ecology classes, this
research demonstrates an injective mapping of VERA’s simulation output onto
common curves found in mathematical ecological modeling. While VERA is a
platform primarily focused on conceptual modeling, making no claim to the
accuracy of its simulations to real world ecology, it is still pedagogically
useful as its outputs are likely to align with the mathematical models
students are learning in their classrooms. VERA is already partially justified
in terms of content validity due to its use in ecology courses and its
relation to the Smithsonian sponsored ecology website Encyclopedia of Life
(EOL). Despite that, this research developed a methodology by which modeling
and simulation based pedagogical tools can be evaluated for content validity,
and using VERA as an example, we were able to show successful overlap between
these methods in larger datasets.

###  5.4 Future Work

The clearest direction for future work is the application of this research to
a larger dataset. This dataset will likely be the entirety of the VERA
database, made of up over 7000 models. A dataset of this size would be the
best foundation for the model comparison part of VERA’s personalized coaching
system.

A proper analysis of outlier curves and their underlying causes in the VERA
system is another path for future work. Patterns in these curves may reveal
either more complex known ecological patterns currently unknown to the VERA
team or it may reveal something corresponding to the software of VERA itself.
Further research into outlier patterns is the best way to differentiate
between these two possible causes, and both lead to a greater understanding of
the VERA system as a whole.

##  6 Conclusion

In this work, we explore analytics as a method by which we can assess content
validity for conceptual modeling and simulation–based tools in domains that
heavily rely on mathematical modeling. As an example we apply the methodology
to user generated content from VERA—an ecology education tool. Using common
methods taken from both analysis of real world data, in the form of
hierarchical clustering, and from mathematical ecological models, in the form
of curve fitting, we analyze the outputs of an agent–based simulation system.
We found the greatest success using the largest data set reaching a 89.38%
success rate when compared to the curve fitting algorithm, with some
indication that an even larger data set would perform better on average. We
also saw a convergence of the initial clustering algorithm success rate, and
the KNN validation step, with the difference between the two evaluation
methods approaching as little as 3.86%. The high rate of success in the data
superset serves as a promising indication that this methodology could be
applied to any other simulation–based education tools for which their domain
has a strong modeling component, such as economics or other social sciences,
and its application would give similar results for content valid systems.

## References

  * [20222] Modeling Ecosystem Dynamics, jun 8 2022.  [Online; accessed 2023-07-16].
  * [ABH+18] Sungeun An, Robert Bates, Jennifer Hammock, Spencer Rugaber, and Ashok Goel.  Vera: popularizing science through ai.  In Artificial Intelligence in Education: 19th International Conference, AIED 2018, London, UK, June 27–30, 2018, Proceedings, Part II 19, pages 31–35. Springer, 2018.
  * [AGG+10] M Anand, A Gonzale, F Guichard, J Kolasa, and L Parrott.  Ecological systems as complex systems: Challenges for an emerging science.  Diversity, 2:395–410, 2010.
  * [AMST11] Wolfgang Aigner, Silvia Miksch, Heidrun Schumann, and Christian Tominski.  Visualization of time-oriented data, volume 4.  Springer, 2011.
  * [BA11] Seely AJ. Bravi A, Longtin A.  Review and classification of variability analysis techniques with clinical applications.  Biomedical engineering online, pages 1–27, 2011.
  * [BSM+21] K. Bogomolova, A. H. Sam, A. T. Misky, C. M. Gupte, P. H. Strutton, T. J. Hurkxkens, and B. P. Hierck.  Development of a virtual three‐dimensional assessment scenario for anatomical education.  Anatomical sciences education, 14, 2021.
  * [CCST15] Ruyin Cao, Jin Chen, Miaogen Shen, and Yanhong Tang.  An improved logistic method for detecting spring vegetation phenology in grasslands from modis evi time-series data.  Agricultural and Forest Meteorology, 200:9–20, 2015.
  * [CH16] D. A. Cook and R. Hatala.  Validation of educational assessments: a primer for simulation and beyond.  Advances in simulation, 1:1–12, 2016.
  * [CZH+14] D. A. Cook, B. Zendejas, S. J. Hamstra, R. Hatala, and R. Brydges.  What counts as validity evidence? examples and prevalence in a systematic review of simulation-based assessment.  Advances in Health Sciences Education, 19, 2014.
  * [DB15] Crystal Day-Black.  Gamification: An innovative teaching-learning strategy for the digital nursing students in a community health nursing course.  ABNF Journal, 26, 2015.
  * [DlTOC+21] R. De la Torre, B. S. Onggo, C. G. Corlu, Nogal M., and A. A. Juan.  The role of simulation and serious games in teaching concepts on circular economy and sustainable energy.  Energies, 14:1138, 2021.
  * [ENB12] Matthew R Evans, Ken J Norris, and Tim G Benton.  Predictive ecology: systems approaches, 2012.
  * [Gam98] S. Gamito.  Growth models and their use in ecological modelling: an application to a fish population.  Ecological modelling, 113:83–94, 1998.
  * [GJ15] A Goel and D Joyner.  Impact of a creativity support tool on student learning about scientific discovery processes.  In Proceedings of the Sixth International Conference on Computational Creativity, 2015.
  * [Got95a] N. J. Gotelli.  A primer of Ecology.  Sinauer Associates Incorporate, 1995.
  * [Got95b] N. J. Gotelli.  A primer of Ecology.  Sinauer Associates Incorporate, 1995.
  * [GRV09] A. K. Goel, S. Rugaber, and S. Vattam.  Structure, behavior, and function of complex systems: The structure, behavior, and function modeling language.  Ai Edam, 23, 2009.
  * [HCR06] M Holcombe, S Coakley, and Smallwood R.  A general framework for agent-based modelling of complex systems.  Proceedings of the 2006 European conference on complex systems, 1, 2006.
  * [J.04] Vandermeer J.  Coupled oscillations in food webs: Balancing competition and mutualism in simple ecological models.  The American Naturalist, 163:857–867, 2004.
  * [JE02] Per Jonsson and Lars Eklundh.  Seasonality extraction by function fitting to time-series of satellite sensor data.  IEEE transactions on Geoscience and Remote Sensing, 40(8):1824–1832, 2002.
  * [JGP14] David A Joyner, Ashok K Goel, and Nicolas M Papin.  Mila–s: generation of agent-based simulations from conceptual models of complex systems.  In Proceedings of the 19th international conference on intelligent user interfaces, pages 289–298, 2014.
  * [JHC74] Gauch Jr, G. Hugh, and Gene B. Chase.  Fitting the gaussian curve to ecological data.  Ecology, 55, 1974.
  * [JMA14] & Waring T. M. Janssen M. A., Lee A.  Experimental platforms for behavioral experiments on social-ecological systems.  Ecology and Society, 19, 2014.
  * [JRKB19] M. Jaxa-Rozen, J. H. Kwakkel, and M. Bloemendal.  A coupled simulation architecture for agent-based/geohydrological modelling with netlogo and modflow.  Environmental modelling & software, 115, 2019.
  * [Kas17] A. Kassambara.  Practical guide to cluster analysis in R: Unsupervised machine learning.  2017\.
  * [KTGA13] Dimitrios Kotsakos, Goce Trajcevski, Dimitrios Gunopulos, and Charu C Aggarwal.  Time-series data clustering., 2013.
  * [KW21] S. C. Kong and Y. Q. Wang.  Item response analysis of computational thinking practices: Test characteristics and students’ learning abilities in visual programming contexts.  Computers in Human Behavior, 122, 2021.
  * [LLE08] D. P. Demaster L. L. Eberhardt, J. M. Breiwick.  Analyzing population growth curves.  Oikos, 117:1240–1246, 2008.
  * [Łuc16] Maciej Łuczak.  Hierarchical clustering of time series data with parametric derivative dynamic time warping.  Expert Systems with Applications, 62:116–130, 2016.
  * [LVVC11] S. Lhermitte, J. Verbesselt, W.W. Verstraeten, and P. Coppin.  A comparison of time series similarity measures for classification and change detection of ecosystem dynamics.  Remote sensing of environment, 115:3129–3152, 2011.
  * [MC17] F. Murtagh and P. Contreras.  Algorithms for hierarchical clustering: an overview.  Wiley Interdisciplinary Reviews: Data Mining and Knowledge Discovery, 2, 2017.
  * [Mes95] S. Messick.  Validity of psychological assessment: Validation of inferences from persons’ responses and performances as scientific inquiry into score meaning.  American psychologist, 50, 1995.
  * [MTD+18] J. L. McGrath, J. M. Taekman, P. Dev, D. R. Danforth, D. Mohan, N Kman, and K. Won.  Using virtual reality simulation environments to assess competence for emergency medicine learners.  Academic Emergency Medicine, 25, 2018.
  * [OES07] Mahamed GH Omran, Andries P Engelbrecht, and Ayed Salman.  An overview of clustering methods.  Intelligent Data Analysis, 11(6):583–605, 2007.
  * [OR06] P. Ormerod and B. Rosewell.  Validation and verification of agent-based models in the social sciences.  International workshop on epistemological aspects of computer simulation in the social sciences, 2006.
  * [oTDL23] Georgia Institute of Technology Design and Intelligence Lab.  Virtual ecological research assistant (vera), 2023.
  * [Ott03] Julio M. Ottino.  Complex systems.  American Institute of Chemical Engineers, 49, 2003.
  * [Pie81] E.C Pielou.  The usefulness of ecological models: a stock-taking.  The Quarterly Review of Biology, 56:17–31, 1981.
  * [PWL+14] Cynthia S Parr, Mr Nathan Wilson, Mr Patrick Leary, Katja S Schulz, Ms Kristen Lans, Ms Lisa Walley, Jennifer A Hammock, Mr Anthony Goddard, Mr Jeremy Rice, Mr Marie Studer, et al.  The encyclopedia of life v2: providing global access to knowledge about life on earth.  Biodiversity data journal, (2), 2014.
  * [S.] An S.  An, sung disseration.  Sung dissertation.
  * [SB13] SB Snider and JN Brimlow.  An introduction to population growth.  Nature Education Knowledge, 4(4):3, 2013.
  * [TKG12] J. C. Thiele, W. Kurth, and V. Grimm.  Rnetlogo: An r package for running and exploring individual‐based models implemented in netlogo.  Methods in Ecology and Evolution, 3, 2012.
  * [TW04a] S. Tisue and U Wilensky.  Netlogo: A simple environment for modeling complexity.  International conference on complex systems, 21:16–21, 2004.
  * [TW04b] Seth Tisue and Uri Wilensky.  Netlogo: A simple environment for modeling complexity.  In International conference on complex systems, volume 21, pages 16–21. Citeseer, 2004.
  * [WB97] Lu F. Winterhalder B.  A forager‐resource population ecology model and implications for indigenous conservation.  Conservation Biology, 11:1354–1364, 1997.
  * [WF90] Barbara Y White and John R Frederiksen.  Causal model progressions as a foundation for intelligent learning environments.  Artificial intelligence, 42(1):99–157, 1990.
  * [WFM07] P. Windrum, G. Fagiolo, and A. Moneta.  Empirical validation of agent-based models: Alternatives and prospects.  Journal of Artificial Societies and Social Simulation, 10, 2007.
  * [ZWX+20] Linglin Zeng, Brian D Wardlow, Daxiang Xiang, Shun Hu, and Deren Li.  A review of vegetation phenological metrics extraction using time-series, multispectral satellite data.  Remote Sensing of Environment, 237:111511, 2020.