# Optimizing the Resolution of Hydrodynamic Simulations for MCRaT Radiative
Transfer Calculations

Jose Arita-Escalante Southeastern Universities Research Association,
Washington, D.C. 20005, USA  Astrophysics Science Division, NASA Goddard Space
Flight Center,Greenbelt, MD 20771, USA  Center for Research and Exploration in
Space Science and Technology, NASA/GSFC, Greenbelt, Maryland 20771, USA  Tyler
Parsotan Center for Space Science and Technology, University of Maryland
Baltimore County, 1000 Hilltop Circle, Baltimore, MD 21250, USA  Astrophysics
Science Division, NASA Goddard Space Flight Center,Greenbelt, MD 20771, USA
Center for Research and Exploration in Space Science and Technology,
NASA/GSFC, Greenbelt, Maryland 20771, USA  S. Bradley Cenko Astrophysics
Science Division, NASA Goddard Space Flight Center,Greenbelt, MD 20771, USA
Joint Space-Science Institute, University of Maryland, College Park, MD 20742,
USA

###### Abstract

Despite their discovery about half a century ago, the Gamma-ray burst (GRB)
prompt emission mechanism is still not well understood. Theoretical modeling
of the prompt emission has advanced considerably due to new computational
tools and techniques. One such tool is the PLUTO hydrodynamics code, which is
used to numerically simulate GRB outflows. PLUTO uses Adaptive Mesh Refinement
to focus computational efforts on the portion of the grid that contains the
simulated jet. Another tool is the Monte Carlo Radiation Transfer (MCRaT)
code, which predicts electromagnetic signatures of GRBs by conducting photon
scatterings within a jet using PLUTO. The effects of the underlying resolution
of a PLUTO simulation with respect to MCRaT post-processing radiative transfer
results have not yet been quantified. We analyze an analytic spherical outflow
and a hydrodynamically simulated GRB jet with MCRaT at varying spatial and
temporal resolutions and quantify how decreasing both resolutions affect the
resulting mock observations. We find that changing the spatial resolution
changes the hydrodynamic properties of the jet, which directly affect the
MCRaT mock observable peak energies. We also find that decreasing the temporal
resolution artificially decreases the high energy slope of the mock observed
spectrum, which increases both the spectral peak energy and the luminosity. We
show that the effects are additive when both spatial and temporal resolutions
are modified. Our results allow us to understand how decreased hydrodynamic
temporal and spatial resolutions affect the results of post-processing
radiative transfer calculations, allowing for the optimization of hydrodynamic
simulations for radiative transfer codes.

Keywords pending

††software: Python (Python Core Team, 2019), Astropy (Astropy Collaboration et
al., 2013, 2018, 2022), NumPy (Harris et al., 2020), MCRaT (Lazzati, 2016;
Parsotan & Lazzati, 2018; Parsotan et al., 2018; Parsotan & Lazzati, 2021,
2022), ProcessMCRaT (Parsotan, 2021)

##  1 Introduction

A number of different theories have been created to understand the phenomena
of Gamma-ray bursts (GRBs) since their initial discovery in the 1960’s
(Klebesadel et al., 1973). One of the earliest models used to explain GRB
prompt emission was the Synchroton Shock Model (SSM)(Rees & Meszaros, 1994),
which considers radiation generated when shells with different Lorentz factors
collide with each other outside of the photospheric region (Daigne et al.,
2011). The collisions of these shells create perturbations of the magnetic
fields that lead to the excitation of leptons which then emit synchroton
radiation. The SSM can naturally explain GRB properties such as lightcurve
variability and the observed nonthermal spectra. Nevertheless, it fails to
agree with observed correlations of GRBs such as the Amati and Yonetoku
relations (Amati et al., 2002; Yonetoku et al., 2004; Zhang & Yan, 2011).

Another model explaining the prompt emission mechanism is the photospheric
model, which explains the phenomenon by describing thermal radiation that
originates deep within a relativistic jet (Rees & Mészáros, 2005). The
radiation is initially in a part of the jet with a high optical depth, leading
to many interactions between the photons and the matter in the jet. As the jet
expands, it becomes optically thin, allowing photons to leave the jet’s
photosphere and travel to the observer without additional interactions with
the GRB jet. The photospheric model is able to reproduce correlations that the
SSM cannot, but is unable to replicate non-thermal spectral low and high-
energy tails without the consideration of the photospheric region
(Beloborodov, 2010; Pe’er, 2008; Pe’er & Ryde, 2011) and subphotospheric
dissipation events (Chhotray & Lazzati, 2015).

With the aid of computational tools, we have been better able to understand
the physics of GRBs. Previous studies have conducted rigorous radiative
transfer calculations, however; they have assumed that the jet structure has
been simplified into an analytic profile (Ito et al., 2013, 2014; Vurm &
Beloborodov, 2016). In contrast to radiative transfer calculations, other
studies have utilized hydrodynamic (HD) calculations to simulate complex jet
structures, but these only provide information about the matter within the jet
(Lazzati et al., 2009, 2013; López-Cámara et al., 2014), which leads to a lack
of information regarding the evolution of the radiation.

The state of the art method to account for both of these assumptions is to
perform post-processing radiative transfer calculations on a hydrodynamic (HD)
simulated jet using Monte Carlo methods. There have been tools developed to
perform post-processing radiative transfer calculations such as the ones
developed by Ito et al. (2015, 2019) and the Monte Carlo Radiation Transfer
(MCRaT) code (Lazzati, 2016; Parsotan & Lazzati, 2018; Parsotan et al., 2018;
Parsotan & Lazzati, 2021, 2022). MCRaT was developed to conduct radiative
transfer calculations on HD simulations to generate mock observations of
simulated GRBs using the photospheric model. The impact that HD simulation
resolutions have on MCRaT post-processing radiative transfer calculations has
not been studied yet. Ensuring that radiative transfer calculations are
converged and accurate is critical to testing GRB prompt emission theories
against observations.

Here, we present an analysis of HD resolution and its effect on post-
processing radiative transfer calculations for simulated GRB mock observables.
Section 2 outlines the code used to create the HD jet, the code used to
perform the radiative transfer calculations, and the way in which resolutions
are quantified. Sections 3 and 4 show the effect that HD simulation
resolutions have on radiative transfer calculations and the physical
implications of these results.

##  2 Methods

In this section, we outline the methods to our analysis. In Section 2.1, we
discuss the codes used in our analysis. In Section 2.2, we quantify
convergences in our simulations.

###  2.1 Codes Used

Here, we discuss the codes used in our study. Section 2.1.1 highlights the
tools used to create the numerical HD GRB simulation. In Section 2.1.2, we
discuss the tools used to conduct post-processing radiative transfer
calculations and analyze the results to generate mock observables for
simulated GRBs.

####  2.1.1 PLUTO

Figure 1: Frame properties at the end of the PLUTO 16TI simulation. These
quantities are read at time given by
$t_{sim}=\frac{\text{frame}}{\text{fps}}=527.6$ s. Panel (a) shows the Lorentz
factor profile at this time. Panel (b) shows the mass density profile at this
time. Panel (c) shows the AMR Grid shape and levels at this time. In this last
panel, we can see the computational efforts are dedicated to the higher
refinement area of the frame, where the jet is located. We note that even
though Panel (c) shows 7 levels of refinement, in our analysis we highlight
the 5 highest levels of spatial refinement in any HD frame.

PLUTO is a numerical solver for systems of partial differential equations in
the context of astrophysical fluid dynamics (Mignone et al., 2007). PLUTO uses
CHOMBO Adaptive Mesh Refinement (AMR) to focus computational efforts on the
most relevant parts of the HD simulation (Mignone et al., 2012). Here, we used
the PLUTO hydrodynamics code with AMR to simulate the propagation of a long
gamma-ray burst (LGRB) jet from a 16TI stellar progenitor, taken from Woosley
& Heger (2006). The stellar progenitor profile was interpolated onto the PLUTO
grid using the code’s capability to do such operations. Following the
prescription provided by Lazzati et al. (2013), we inject the jet with a
constant luminosity of $5.33\times 10^{50}$ erg $\text{s}^{-1}$ for 100 s from
an injection radius of $1\times 10^{9}$ cm, with an initial Lorentz factor of
5, an opening angle $\theta_{0}=10^{\circ}$, and an internal over rest-mass
energy ratio, $\eta=80$. The simulation is given in a 2D spherical coordinate
system. The simulation domain in PLUTO is logarithmic in radius from $1\times
10^{9}$ cm to $5.6\times 10^{14}$ cm, with 1600 cells in the radial direction,
although the simulation is only carried out until the jet head reaches $\sim
2\times 10^{13}$ cm. It also extends in polar angle from $0^{\circ}$ to
$90^{\circ}$ covered within 160 cells. The AMR refinement is set such that the
jet is followed with a resolution of at least $1\times 10^{9}$ cm along the
jet axis, leading to a maximum of 7 refinement levels needed for the AMR. The
state of the jet is saved with a frame rate of 5 frames per second.

A snapshot of the simulation can be found in Figure 1. These plots show a
snapshot taken at the end of the simulation, at $t_{sim}=527.6$ s. Figure 1(a)
shows the Lorentz factor as a function of position in the frame. Figure 1(b)
shows the mass density profile of this last frame. Figure 1(c) shows the AMR
grid for this frame where level 1 is the base hydrodynamic grid that was
previously defined. The highest level of refinement is level 7 in this figure
which is focused on the core of the jet.

We show the refinement and the convergence of the hydrodynamic properties of
the simulation at the initial moment of photon injection and at the end of the
simulation in Figure 2.

We select a shell of grid cells around a radius of $1.3\times 10^{12}$ cm at
50 s in the simulation and show the spatial resolution that is achieved at
each refinement level at that time. We also show the convergence of the bulk
Lorentz factor, $\Gamma$, density, $\rho$, and temperature , $T$111Temperature
is calculated assuming that the jet is radiation dominated., as we traverse
the different refinement levels. This convergent behavior is present
throughout the whole simulation. In Figure 2 we also show the HD properties of
a shell of grid cells located at $1.5\times 10^{13}$ cm at 527.6 s in the
simulation, which is the last frame of our simulation. As is outlined in
Section 2.1.2 we focus on the 5 highest spatial refinement levels in the
simulation regardless of the hydrodynamic simulation frame that is read in.
Thus, the final frame quantities for refinement level 5 shown in Figure 2
correspond to the black grid representing level 7 shown in 1(c).

Figure 2: HD properties at different refinement levels for the PLUTO 16TI
simulation as functions of refinement level at the first moment of photon
injection and final frame of the simulation. We took the initial frame
measurements of a shell located in radius $r=1.3\times 10^{12}$ cm and time
$t_{sim}=\frac{\text{frame}}{\text{fps}}=50$ s. We took the final frame
measurements of a shell located in radius $r=1.5\times 10^{13}$ cm and time
$t_{sim}=\frac{\text{frame}}{\text{fps}}=527.6$ s. Panel (a) shows the average
HD cell radius size ($\Delta r$) in cm. Panel (b) shows the average density in
$\text{g cm}^{-3}$. Panel (c) shows the average temperature in K. Panel (d)
shows the average bulk Lorentz factor. The early-time conditions are shown as
red triangle markers and the final conditions are shown in blue circle
markers. Throughout our analysis we focus on the 5 highest spatial refinement
levels in a given HD simulation frame. Thus, the final frame quantities for
refinement level 5 shown in this figure correspond to the black grid
representing level 7 shown in 1(c).

####  2.1.2 MCRaT and ProcessMCRaT

The MCRaT222https://github.com/lazzati-astro/MCRaT code conducts radiative
transfer calculations to compute the electromagnetic (EM) signature of HD
simulated GRB jets. MCRaT reads in HD simulations of GRB jets and performs
Compton scatterings between the injected photons and matter in the jet. The
code first reads in the simulation snapshots and injects photons within a
specified region of the jet, where the number of injected photons depends on
the comoving temperature of the jet. These photons are then propagated between
and within each subsequent snapshot of the HD simulation 333There are no
discretized emission regions for photons in the MCRaT calculations as there
are in other systems where monte carlo radiative transfer methods are used
(e.g. stars in Stamatellos & Whitworth (2005) or Parsotan et al. (2021)).
There are mechanisms that allow for continuous photon emission such as that
from cyclo-synchrotron emission (Parsotan & Lazzati, 2021) however we do not
consider that mechanism in this study. Additionally, subphotospheric photon
mediated shocks are expected to produce photons in GRB jets but these physics
are not yet captured in large scale, global HD or radiative transfer
simulations (Ito et al., 2018). The post-processing radiative transfer method
that MCRaT employs is necessary for conducting photon transport calculations
in special relativity within a medium that is time dependent.

MCRaT can run two different radiative transfer calculations. The first one is
based off of reading in an HD numerical simulation of a GRB jet. In our study,
we use the PLUTO 16TI simulation mentioned in Section 2.1.1. The outflow given
by this numerical simulation introduces numerically induced errors in the HD
properties of the grid, so the effect of HD resolution on the post-processing
radiative transfer calculations can be assessed. The other type of radiative
transfer calculation MCRaT is capable of running is one of an analytic
spherical outflow, or fireball. In the spherical outflow case, MCRaT takes the
HD simulation files and overwrites the HD properties with those of an analytic
outflow with only outward radial velocity components. This analytic outflow is
a function of cell radius in which the spherical outflow is accelerated until
an asymptotic Lorentz factor is reached. By using an analytic spherical
outflow, we can understand how just the HD resolution has an effect on MCRaT
mock observables.

We set up our spherical outflow to have an asymptotic Lorentz factor
$\Gamma_{\infty}=100$, luminosity $L=10^{54}\text{ erg}$ $\text{s}^{-1}$ and
saturation radius $r_{0}=10^{8}\text{ cm}$.

Spatial resolution for the HD simulation takes the form of various AMR
refinement levels. Since the PLUTO AMR simulation that we conducted
dynamically changes the number of refinement levels in order to maintain a
resolution element size of $\sim 10^{9}$ cm, MCRaT reads in the $n$th highest
refinement level at any given frame. Thus, level 5 is the highest refinement
level at any time in the HD simulation, level 4 would be the second highest
level, all the way until level 1, which is the lowest refinement level. We run
multiple MCRaT simulations where we specify which of these spatial refinement
levels should be read in from the PLUTO simulation.

As outlined in Section 2.1.1, the PLUTO simulation we use has a framerate of 5
frames per second (fps). For the context of our analysis, this is our highest
temporal refinement level. We artificially vary the frame rate of our PLUTO
simulation by telling MCRaT to only read every $n$th frame. With this method,
we can achieve our desired framerate while maintaining the total simulation
time for the HD simulation. In order to keep the same simulation time, each
step in time ($\Delta t_{sim}$) has a specific number of frames assigned to
it. As the resolution is lowered, the first and last frame in each $\Delta
t_{sim}$ stay the same while the number of frames in between these two is
lowered. This leads to a “choppy” simulation. In order to keep the simulation
realistic as a function of time, each pair of subsequent frames varies more
significantly as the temporal resolution is lowered. Reducing the framerate by
a factor of 2 each time would be analogous to contiguous spatial refinement
level HD cell radius increasing by a factor of 2. Therefore, 5 fps would be
analogous to spatial refinement level 5, 2.5 fps would be analogous to spatial
refinement level 4 and so on. This gives us a way to align the spatial
refinement levels and come up with a clear way to mix and match temporal and
spatial refinement levels, allowing us to investigate the effects of these
changes combined with and independent of one another.

For all our MCRaT simulations, we kept our parameters as constant as possible.
We injected photons into the HD simulation at an angle range of
$0^{\circ}-9^{\circ}$ and radius at $\sim 10^{12}$ cm. We simulate photons
within the first $100$ s of the PLUTO 16TI simulation, the time for which the
GRB jet is active. Additionally, we simulated $\sim 10^{5}-10^{6}$ photons per
MCRaT simulation.

In order to analyze the output of MCRaT’s simulations, we used
ProcessMCRaT444https://github.com/parsotat/ProcessMCRaT (Parsotan, 2021).
ProcessMCRaT is a Python library developed to analyze and manipulate the
output of MCRaT radiative transfer calculations. ProcessMCRaT fits the mock
observed spectrum with a Band function (Band et al., 1993) to calculate its
low and high energy slopes, $\alpha$ and $\beta$ respectively, and its break
energy, $E_{\text{0}}$. The peak energy of the spectrum can then be calculated
as $E_{\text{pk}}=(2+\alpha)E_{\text{0}}$ . ProcessMCRaT also has the
capability to create mock lightcurves for MCRaT simulated GRBs. ProcessMCRaT
uses the photons in the last HD frame of the simulation to calculate mock
observables, which corresponds to $t_{sim}=527.6$ s for the PLUTO simulation
used in this study which is shown in Figure 1. The photons are assumed to be
in an optically thin regime of the jet and thus are no longer scattering with
the matter in the jet, changing their energies, as they propagate to the
observer. Choosing to construct mock observables at an earlier time in the
simulation would lead to incorrect transient mock observables instead of
steady state mock observables expected from photons that are no longer coupled
to the jet.

In producing our mock observables, we placed a mock observer at
$r_{\text{obs}}=10^{14}$ cm at various angles
$\theta_{\text{obs}}=1^{\circ},3^{\circ},5^{\circ}$ and $8^{\circ}$ from the
GRB jet axis. The opening angle for the area in which the observer detects
photons in $\Delta\theta_{\text{obs}}=4^{\circ}$. We set the spectral fit for
the observables to be that of a Band function including all photons at an
energy range of $0.1$ \- $4000$ keV.

We numerically integrated spectra with respect to energy to get luminosities,
$L_{\mathrm{iso}}$. We also numerically integrated lightcurves with respect to
time to obtain total isotropic energies, $E_{\mathrm{iso}}$.

###  2.2 Quantifying Convergence Within Radiative Transfer Calculations

Since there are two dimensions of change in refinement (spatial and temporal),
we populate a $5\times 5$ matrix that contains entries for its spatial and
temporal refinement levels.

In order to quantify convergence in MCRaT mock observables between one
spatial/temporal level and another, we define the percent change variable
$\zeta_{\text{Prop}}$ as:

$\displaystyle\zeta_{\text{Prop}}(\text{lev}(n),$ $\displaystyle\text{fps})=$
$\displaystyle\left|\frac{\text{Prop}(\text{lev}(n),\text{fps})-\text{Prop}(\text{lev}(5),5\text{
fps})}{\text{Prop}(\text{lev}(5),5\text{ fps})}\right|.$ (1)

Equation 2.2 represents a comparison of any level of temporal and spatial
refinement with the highest combination of these both (spatial refinement 5
and 5 fps in our context) for any particular property (called Prop in Equation
2.2) of the GRB EM signature, such as $\alpha$, $E_{\mathrm{pk}}$, or
$L_{\mathrm{iso}}$555The calculations can be found here:
https://doi.org/10.5281/zenodo.8139769.

This gives us a way to quantify the deviation at each level compared to the
highest level for each GRB EM property. For the analysis of our results, the
quantity $\zeta_{\text{Prop}}$ will be used to quantify deviations in
resulting mock observables at different refinement levels.

##  3 Results

Here, we outline the results of our findings for the spherical outflow case
and the 16TI simulation as described in Section 2. Our analysis shows the same
trends for our mock observer angle $\theta_{\text{obs}}$ at all angles
mentioned in Section 2.1.2. As a result, this section will only focus on
$\theta_{\text{obs}}=1^{\circ}$.

###  3.1 Spherical Outflow

####  3.1.1 Spectra

Figure 3: Spectra of a spherical outflow for different refinement levels. The
solid purple line is a blackbody spectrum that peaks at each spectrum set’s
highest refinement level. Panel (a) shows the spectra of a spherical outflow
profile, but at different spatial refinement levels while maintaining the same
highest temporal resolution constant. Panel (b) shows the spectra generated
with analytic outflow at different temporal resolutions, maintaining the
highest spatial resolution constant. Panel (c) shows the spectra with matching
temporal and spatial resolution levels.

Figure 3 shows spectra of a spherical outflow at different spatial and
temporal resolutions. Figure 3(a) shows spectra at the highest temporal
resolution and varying spatial resolutions. Figure 3(b) shows spectra at the
highest spatial resolution and varying temporal resolutions. Figure 3(c) shows
spectra at matching temporal and spatial resolution levels.

When reducing the spatial resolution, we see an artificial increase of the
peak energy of the spectrum. The higher HD cell sizes in lower resolutions
cause the injection coordinates for photons to have different HD properties.
As seen in Figure 2, lower spatial resolution levels have higher temperatures,
which means the injected photons will have higher energies, causing a higher
$E_{\mathrm{pk}}$. Figure 4 shows the spectral peak energies at varying
temporal and spatial resolutions. The effect of lowering the spatial
resolution on spectral peak energies can be seen in Figure 4(a).

Figure 4: Spectral Peak energies of a spherical outflow for different
refinement levels. Panel (a) shows the $E_{\mathrm{pk}}$ of the same spherical
outflow simulation, but at different spatial refinement levels while
maintaining the same highest temporal resolution constant. Panel (b) shows the
$E_{\mathrm{pk}}$ generated with the same analytic outflow simulation at
different temporal resolutions, maintaining the highest spatial resolution
constant. Panel (c) shows the $E_{\mathrm{pk}}$ with matching temporal and
spatial resolution levels. The error bars in panels (a), (b) and (c) are
present, but are encompassed within the markers.

Reducing the spatial resolution causes an increase in the luminosity of the
lightcurves. This happens since there is now more energy in the spectrum. The
spectra then are shifted up in luminosity and to the right in energies, while
still maintaining a blackbody shape. This effect can also be seen in Figure
4(a).

Reducing the temporal resolution of the simulation does not affect the
spectrum of the spherical outflow in any significant manner. Since there is no
change in spatial resolution, the injected photons read in the same HD values
regardless of the temporal resolution. This effect can be seen in Figure 3(b).
For this reason, the peak energies are not affected as the temporal resolution
is decreased, as seen in Figure 4(b). The luminosity of the spectra seems to
be slightly decreased as the temporal resolution decreases. There is not a lot
of variation in the spectral shape and properties since the homologous
expansion present in a spherical outflow does not depend on time.

Mixing spatial and temporal resolutions shows similar trends to only changing
spatial resolutions. This is to be expected, since the analytic spherical
outflow is defined to be time-independent. This effect can be seen in Figure
3(c).

Analyzing the peak energies $E_{\mathrm{pk}}$ for different refinement levels
confirms what we observed in the spectra in Figure 3. As seen in Figure 4(a),
there is an increase in peak energy as the spatial refinement level is
decreased. This is seen as a shift to the right in the spectra in Figure 3(a).
For temporal resolutions, as seen in Figure 3(b), there is no significant
change between levels. This can be seen in the $E_{\mathrm{pk}}$ values in
Figure 4(b) where the values are similar to one another. When mixing spatial
and temporal resolutions, there is an additive behavior in the differences at
various levels.

One other important quantity to observe is the luminosity at different spatial
and temporal refinement levels. Figure 5 shows the spectral luminosity at
varying spatial and temporal resolutions with a spherical outflow. For
differing spatial refinement levels, there is an artificial increase in
luminosity as seen as a shift upwards in the spectra in Figure 3(a). This
effect can be seen in Figure 5(a). For differing temporal resolutions,
luminosities tend to oscillate not too far way from each other as seen in
Figure 5(b). This aligns with the very similar spectra seen in Figure 3(b).
When there is a mix of temporal and spatial refinement levels, luminosities
also are artificially increased as the resolution goes down. This is the
effect of the same phenomenon happening in spatial resolutions. The additive
effect of combining spatial and temporal resolutions can be seen in Figure
5(c).

Figure 5: Spectral luminosities of a spherical outflow for different
refinement levels. Panel (a) shows the luminosities using the same analytic
outflow simulation, but at different spatial refinement levels while
maintaining the same highest temporal resolution constant. Panel (b) shows the
luminosities generated with the same spherical outflow simulation at different
temporal resolutions while holding the highest spatial resolution constant.
Panel (c) shows the luminosities with matching temporal and spatial resolution
levels. The error bars in panels (a), (b) and (c) are present, but are
encompassed within the markers.

####  3.1.2 Lightcurves

Figure 6 shows lightcurves of a spherical outflow at different spatial and
temporal resolutions. Figure 6(a) shows lightcurves at the highest temporal
resolution and varying spatial resolutions. Figure 6(b) shows lightcurves at
the highest spatial resolution and varying temporal resolutions. Figure 6(c)
shows lightcurves at matching temporal and spatial resolution levels.

In a spherical outflow, the lightcurves at different temporal and spatial
resolutions have roughly the same shape. If the spatial resolution is
decreased, as seen in Section 3.1.1, the luminosity of the spectrum is
increased by an upwards shift of the spectrum. This causes lightcurves to have
a higher luminosity. Analytic spherical outflows should have relatively
“constant” lightcurves. This can be seen in Figure 6(a).

Reducing the temporal resolution does not change the luminosity of the
lightcurves like spatial resolutions do. All temporal resolutions seem to
roughly have the same luminosity. Analyzing the lightcurves, there is an
increase in variability as the lightcurves oscillate around one “average”
value of the lightcurve. This variation is due to the fact that the photons
are not being smoothly injected in a thin shell. Instead, the MCRaT algorithm
has to determine which HD cells are the most energetic within a larger set of
HD cells and correspondingly place more photons in those photon dense regions
of the HD simulation. This leads to us only probing the portions of the
outflow with the largest energies. This causes us to no longer get a smooth
stream of photons that are detected as a function of time. This effect can be
seen in Figure 6(b).

Like with the spectra, varying both resolutions at the same time has an
additive effect on the changed properties of the lightcurves. Figure 6(c)
shows how decreasing both the temporal and spatial resolution increases the
luminosity of the lightcurve as well as the variability in the form of an
oscillation around the “average” value of the lightcurve.

Figure 6: Lightcurves of simulated spherical outflow for different refinement
levels. Panel (a) shows the lightcurves using the same analytic outflow
simulation, but at different spatial refinement levels while maintaining the
same highest temporal resolution constant. Panel (b) shows the lightcurves
generated with the same spherical outflow simulation at different temporal
resolutions, maintaining the highest spatial resolution constant. Panel (c)
shows the lightcurves with matching temporal and spatial resolution levels.

###  3.2 16TI HD Simulation

####  3.2.1 Spectra

Figure 7: Spectra of a 16TI HD simulated GRB at different refinement levels.
The solid purple line is a blackbody spectrum that peaks at each spectrum
set’s highest refinement level. Panel (a) shows the spectra of a simulated GRB
using the same PLUTO 16TI simulation, but at different spatial refinement
levels while maintaining the same highest temporal resolution constant. Panel
(b) shows the spectra generated with the same 16TI simulation at different
temporal resolutions, maintaining the highest spatial resolution constant.
Panel (c) shows the spectra of a simulated GRB with matching temporal and
spatial resolution levels.

Figure 7 shows spectra of a GRB simulated with a 16TI stellar progenitor model
at different spatial and temporal resolutions. Figure 7(a) shows the GRB
spectra at the highest temporal resolution and varying spatial resolutions.
Figure 7(b) shows GRB spectra at the highest spatial resolution and varying
temporal resolutions. Figure 7(c) shows GRB spectra at matching temporal and
spatial resolution levels.

Changing spatial resolutions artificially increases the spectral luminosities
and peak energies at the highest framerate, 5 fps, as is shown in Figure 7(a).
We see a less pronounced artificial increase in these quantities simulations
where the framerates is held constant at a lower value then 5 fps. The
introduction of a lower spatial resolution leads to a sudden change in the HD
properties of each grid cell with respect to its neighbors. This abrupt change
in HD properties makes it hard for the photons to maintain equilibrium with
the medium and leads to them being artificially upscattered. Decreased spatial
resolutions lead to higher temperatures which causes an upwards shift of the
spectrum, while keeping the spectral shape relatively constant. This effect is
decreased as the spatial resolution is increased since the smoother HD
behavior makes it easier for the photons to stay in equilibrium with the jet.
The artificially increased spectral high-energy tail due to lower spatial
resolutions and spurious upscatterings can be seen in Figure 7(a).

Lower temporal resolutions also lead to an abrupt change in the HD properties
of the simulation. The abrupt change due to lower temporal resolutions is
different in nature than the change in spatial resolutions, although the end
result is the same. The lower framerate leads to photons scattering in the
same HD frame for longer periods of time. Once the next frame is reached, the
gradient in the HD properties is more pronounced, leading to the photons
upscattering to higher energies. The artificially increased spectral high-
energy tail due to lower temporal resolutions can be seen in Figure 7(b).

Each spectrum showed in Figure 7(b) is at the highest spatial resolution
possible which defines the temperature of the outflow and the number of
photons that are injected into the simulation. Thus, the total number of
photons must be conserved in the spectrum. The dip seen in the spectra of
Figure 7(b) as we lower the temporal resolution is due to the photons near the
peak of the spectrum, at $\sim 100$ keV, being upscattered into the high
energy tail of the spectra. At high temporal refinement levels, this dip in
the spectrum at $\sim 100$ keV does not exist since the photons maintain
thermal equilibrium with the jet.

Combining changes in both temporal and spatial resolutions leads to this
effect being additive. There is upscattering due to large gradients in the
jet’s properties in both space and time however, we find that time resolution
is the dominant factor in accounting for deviations from higher resolutions.
This effect can be seen in Figure 7(c).

####  3.2.2 Lightcurves

Figure 8 shows lightcurves of a GRB simulated with a 16TI stellar progenitor
model at different spatial and temporal resolutions. Only the first $\sim 10$
seconds of the lightcurve are shown to emphasize the effects of lowering
temporal and/or spatial resolutions. Figure 8(a) shows the GRB lightcurves at
the highest temporal resolution and varying spatial resolutions. Figure 8(b)
shows GRB lightcurves at the highest spatial resolution and varying temporal
resolutions. Figure 8(c) shows the GRB lightcurves at matching temporal and
spatial resolution levels.

Decreasing the spatial resolution for the outflow given by the 16TI simulation
has a similar effect to that of decreasing the spatial resolution for a
spherical outflow simulation. Because of the increase in the high energy tail
of the spectrum, there is an increased luminosity that becomes more pronounced
as the spatial resolution is decreased. This can be seen in Figure 8(a).

Changing the temporal resolution for 16TI HD simulated GRBs affects the
lightcurve differently than doing so in a spherical outflow. This is due to
the time dependence now present in the simulated GRB jet. Since the jet
changes more drastically from frame to frame, there is an artificially
enhanced high-energy tail in the spectrum, leading to a higher luminosity in
the lightcurve. Not only is there a higher luminosity but the variation
present in lower temporal resolutions is also present. This variation can be
seen in Figure 8(b).

Changing both temporal and spatial resolutions leads to the combination of
their individual effects. There is an increase in luminosity and in the
variability of the lightcurve. The presence of the increasingly high
lightcurve luminosity and variability can be seen in Figure 8(c).

Figure 8: Lightcurves of the 16TI HD simulated GRB for different refinement
levels. Only the first $\sim 10$ seconds of the simulation are shown to better
visualize the qualitative effects of reducing spatial and/or temporal
resolutions. Panel (a) shows the lightcurves of a simulated GRB using the same
PLUTO 16TI simulation, but at different spatial refinement levels while
maintaining the same highest temporal resolution constant. Panel (b) shows the
lightcurves generated with the same 16TI simulation at different temporal
resolutions, maintaining the highest spatial resolution constant. Panel (c)
shows the lightcurves of a simulated GRB with matching temporal and spatial
resolution levels.

####  3.2.3 Errors in Mock Observables

Figure 9: Matrices containing $\zeta_{\text{Prop}}$ for different EM
properties of the 16TI simulated GRB. Panel (a) shows the deviation
$\zeta_{L_{\text{iso}}}$ for all spatial and temporal refinement levels. Panel
(b) shows the deviation $\zeta_{E_{\text{pk}}}$. Panel (c) shows the deviation
$\zeta_{\alpha}$. Panel (d) shows the deviation $\zeta_{\beta}$.

To better visualize the two degrees of freedom for the change in resolution,
we use a color gradient matrix for all combinations of spatial and temporal
resolutions. Darker colors show larger deviations from the highest level of
refinement, both spatial and temporal. We find the biggest change in the
lowest spatial and temporal resolutions. As we get closer to the highest
resolution (refinement level 5, 5 fps), we find a smaller change in the error.

Luminosity shows the smallest deviation from the highest resolution at $17\%$,
in the (level 4, 5 fps) entry in the resolution matrix. As we deviate from the
highest resolution value in the lower right corner, we start seeing our
deviation change more drastically. The largest deviation comes from the (level
1, 0.3125 fps) entry in the matrix, at a $\sim 5000\%$ change. This result is
expected, since this is the lowest resolution case in both space and time.

In Figure 9(a), we see $\zeta_{L_{\text{iso}}}$ for all refinement levels in
the resolution matrix. The isotropic luminosity comes from integrating the
spectrum with respect to energy. By looking at Figure 7, we can see a clear
reason why these values are changing the way that they are. The high energy
tails present at lower spatial and temporal resolutions give this integral a
higher value. Since there is a larger change in the spectral shape with
respect to temporal resolutions, we see more similar values between spatial
resolutions at constant framerates. We see an identical trend with isotropic
energies $E_{\text{iso}}$. This is expected since this is an integration of
the lightcurve with respect to time.

Peak energies see deviations at around $6-55\%$, as long as we stay within the
lower right $3\times 3$ block in our matrix. As we step outside of this block
(lower than level 3 or 1.25 fps), we start seeing deviations of up to $83\%$.

Figure 9(b) shows $\zeta_{E_{\text{pk}}}$ for all refinement levels in the
resolution matrix for 16TI HD simulated GRBs. There is a present trend in
which there is more variation present as we get further away from the lower
right (level 5, 5 fps) entry. This implies that both spatial and temporal
resolutions affect the resulting peak energies $E_{\text{pk}}$. The spatial
resolution effect is due to the change in temperatures between levels seen in
Figure 2. Since Equation 2.2 calculates the magnitude of the deviation and not
the direction of them, these figures do not encompass this information. The
peak energy values for different refinement levels are decreased as the
spatial resolution is decreased.

When lowering spatial and temporal resolutions, the low energy slope, $\alpha$
varies as low as $0.05\%$ and as high as $67\%$. If we stay in the lower right
$3\times 3$ block of the matrix, we see a deviation of around $1-13\%$. As we
step outside of this block, (lower than level 3 or 1.25 fps), we start seeing
higher deviations. The largest deviation comes from the (level 4, 0.3125 fps)
entry in the matrix, at a $67\%$ change.

Figure 9(c) shows $\zeta_{\alpha}$ for all refinement levels in the resolution
matrix. Figure 7 shows some of the spectra that were evaluated here, but the
trend of change in the low-energy slope $\alpha$ is not as clear by just
looking at the plotted spectra. This figure shows that the lower energy slope
is conserved as resolutions change, showing that lower spatial and temporal
resolutions do not have a large impact on the lower energy subset of photons.
Due to the scatter seen at low energies in the spectra from the low temporal
resolution simulations at 0.3125 fps, as shown in Figure 7(b) where there are
large steps in the spectrum $\lesssim 2$ keV, the fitted $\alpha$ values are
not well constrained which leads to the outliers particularly at the (level 4,
0.3125 fps) entry in the matrix. Additionally, since $\alpha$ and
$E_{\text{pk}}$ are related, the low temporal resolution simulations also
produce outliers in Figure 9(b) see for example the (level 2, 0.3125 fps)
entry of that matrix.

When lowering spatial and temporal resolutions, the high energy slope, $\beta$
is impacted more than the low energy slope $\alpha$. From Figure 7, we can see
that the temporal resolution affects $\beta$ more than the spatial resolution.

Figure 9(d) shows $\zeta_{\beta}$ for all refinement levels in the resolution
matrix. There is an decrease in the high energy slope $\beta$ as the temporal
resolutions are decreased. $\beta$ shows similar values between spatial
resolutions when holding temporal resolutions constant. This effect is due to
photons being upscattered to higher energies as they are shocked by
drastically changing jet properties as new HD frames are loaded. The fitted
Band spectrum attempts to account for the higher amount of high-energy photons
by making the high energy slopes flatter.

##  4 Discussions

We have used MCRaT and PLUTO to quantify the effect that a simulation’s
spatial and temporal resolution has on post-processing radiative transfer
calculations.

We show that the lower spatial and temporal resolutions affect both the shape
of the GRB’s spectrum and lightcurve. The presence of high-energy spectral
tails lead to an increase in isotropic luminosities as resolutions are
decreased. Additionally, lower spatial resolutions lead to higher
temperatures, causing photon energies to be increased, leading an increase in
the normalization of the mock observed spectra.

With the lower spatial and temporal resolutions, photons propagate and scatter
within a choppier HD simulation. As seen in Figure 2, lower spatial
resolutions have higher temperatures, causing injected photons to have
slightly higher energies and increased normalization for their spectrum.

As temporal resolution is decreased, the lightcurves face an increased
variability, and the smoothness of this mock observed property is gone. This
is due to the fact that lower framerates lead to photons not being injected as
smoothly as a function of time. Higher framerates lead to a more continuous
injected photon flow, while lower framerates have more spaced out photon
injections. This leads to the lightcurve having more peaks and troughs
compared to the higher temporal resolution simulations.

Figure 10: $\frac{1/\Gamma^{2}}{\Delta r/r}$ at the final frame of the 16TI HD
simulation. Quantities $>>1$ mean that detected photons are sufficiently able
to probe the smallest spatial scales of the jet. Figure 11: Lightcrossing over
HD cell size ratio at beginning of the simulation. Panel(a) shows
$\frac{c/fps}{\Delta r}=\frac{c\Delta t}{\Delta r}$ at the first moment of
photon injection ($t_{sim}=50$s). Panel (b) shows $\frac{c/fps}{\Delta
r}=\frac{c\Delta t}{\Delta r}$ at the end of the simulation
($t_{sim}=527.6$s).

At the end of the HD simulation, photons have ideally reached the photosphere.
At this point, photons are expected to scatter into an angle of
$\theta=1/\Gamma^{2}$ (Lazzati, 2016). Figure 2(a) shows the radial resolution
at the photon position $r=1.5\times 10^{13}$cm, and Figure 2(d) Shows the bulk
Lorentz factor at this radius. The quantity $\Delta
r/r=\sin{\theta}\approx\theta$ can be used to compare to
$\theta=1/\Gamma^{2}$. If the HD resolution is too low, the HD cell sizes are
larger, which means that the mock observed properties of the photons are
inaccurate due to the detected photons not being able to fully probe the
angular scale visible to the observer. Instead photons are only probing the
properties of just one large HD cell. Figure 10 shows how at higher refinement
levels, the ratio $\frac{1/\Gamma^{2}}{\Delta r/r}$ is large enough to easily
probe the visible portion of the GRB jet regime since this angle encompasses
multiple HD cells. At lower resolutions, we get closer but larger than unity,
as seen for levels 2 and 3. Once level 1 is reached, our ratio is smaller than
1. At this resolution, the photons are not properly probing the visible
angular jet region as they only probe one HD cell.

If we consider the lightcrossing defined by $c/$fps = $c\Delta t$, we get the
distance traveled by a photon in between two frames at any given framerate.
When comparing this distance to the size of HD cells at any given time,
optimally we want the photons to scatter within multiple HD elements as they
move within a given frame. For this to be true, the ratio $\frac{c\Delta
t}{\Delta r}$ must be larger than unity. This would mean that the photons are
able to travel through more than one HD cell during the duration of any given
frame. This allows the photons to properly probe the GRB jet properties.

Figure 11(a) shows the ratio $\frac{c\Delta t}{\Delta r}$ at the beginning of
the 16TI simulation. We can see that higher resolutions with lower framerates
lead to photons probing multiple cells in each individual frame as they travel
through the jet. The lower left corner of the matrix shows framerates and
resolutions in which only one cell is probed in each frame. This phenomenon is
also present and amplified at the end of the simulation, when the HD cell size
is increased from the presence of higher refinement levels. This can be seen
in Figure 11(b).

We can make an analogy with polar coordinates where the photon propagates both
radially, as it diffuses outward with the jet, and in polar angle, where the
observer sees some amount of photons that are propagating towards them if the
observer is located within the photons’ local $1/\Gamma^{2}$ angle. As the
MCRaT photon travels outward within a single simulation frame time step, we
want it to interact with many cells of the GRB jet. This will allow the photon
to change its properties as the jet’s HD properties change radially; see for
example the non-thermal spectra obtained by Parsotan et al. (2018). Once at
the photosphere, we want to receive many photons that are properly probing the
angular size of the jet that the observer is able to see. This can lead to the
appearance of non-thermal spectra, such as the multicolor blackbody (Pe’er &
Ryde, 2011), which is typically seen in GRBs. This can only be seen if the
photons are able to interact with and probe the properties of many different
HD fluid elements in the polar angle direction.

HD simulation resolution needs to be accounted for to retrieve accurate post-
processing radiative transfer calculations. When deciding the HD simulation
resolution for these calculations, one must look at multiple factors to make
an informed decision. These factors include limitations such as storage and
computational resources. Another factor that needs to be considered is what
margin of error is acceptable for one’s particular analysis. In order to save
time, storage and computational resources, if a particular analysis allows for
it, a lower resolution can be chosen, assuming the loss of accuracy in the
mock EM observables is acceptable.

We thank the scientific and data editors for their helpful comments which have
improved the presentation in this paper. The material is based upon work
supported by NASA under award number 80GSFC21M0002. Resources supporting this
work were provided by the NASA High-End Computing (HEC) Program through the
NASA Advanced Supercomputing (NAS) Division at Ames Research Center.

## References

  * Amati et al. (2002) Amati, L., Frontera, F., Tavani, M., et al. 2002, A&A, 390, 81, doi: 10.1051/0004-6361:20020722
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33, doi: 10.1051/0004-6361/201322068
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123, doi: 10.3847/1538-3881/aabc4f
  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, apj, 935, 167, doi: 10.3847/1538-4357/ac7c74
  * Band et al. (1993) Band, D., Matteson, J., Ford, L., et al. 1993, ApJ, 413, 281, doi: 10.1086/172995
  * Beloborodov (2010) Beloborodov, A. M. 2010, MNRAS, 407, 1033, doi: 10.1111/j.1365-2966.2010.16770.x
  * Chhotray & Lazzati (2015) Chhotray, A., & Lazzati, D. 2015, ApJ, 802, 132, doi: 10.1088/0004-637X/802/2/132
  * Daigne et al. (2011) Daigne, F., Bošnjak, Ž., & Dubus, G. 2011, A&A, 526, A110, doi: 10.1051/0004-6361/201015457
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357, doi: 10.1038/s41586-020-2649-2
  * Ito et al. (2018) Ito, H., Levinson, A., Stern, B. E., & Nagataki, S. 2018, Monthly Notices of the Royal Astronomical Society, 474, 2828
  * Ito et al. (2015) Ito, H., Matsumoto, J., Nagataki, S., Warren, D. C., & Barkov, M. V. 2015, ApJ, 814, L29, doi: 10.1088/2041-8205/814/2/L29
  * Ito et al. (2019) Ito, H., Matsumoto, J., Nagataki, S., et al. 2019, Nature Communications, 10, 1504, doi: 10.1038/s41467-019-09281-z
  * Ito et al. (2014) Ito, H., Nagataki, S., Matsumoto, J., et al. 2014, ApJ, 789, 159, doi: 10.1088/0004-637X/789/2/159
  * Ito et al. (2013) Ito, H., Nagataki, S., Ono, M., et al. 2013, ApJ, 777, 62, doi: 10.1088/0004-637X/777/1/62
  * Klebesadel et al. (1973) Klebesadel, R. W., Strong, I. B., & Olson, R. A. 1973, ApJ, 182, L85, doi: 10.1086/181225
  * Lazzati (2016) Lazzati, D. 2016, ApJ, 829, 76, doi: 10.3847/0004-637X/829/2/76
  * Lazzati et al. (2009) Lazzati, D., Morsony, B. J., & Begelman, M. C. 2009, ApJ, 700, L47, doi: 10.1088/0004-637X/700/1/L47
  * Lazzati et al. (2013) Lazzati, D., Morsony, B. J., Margutti, R., & Begelman, M. C. 2013, ApJ, 765, 103, doi: 10.1088/0004-637X/765/2/103
  * Lazzati et al. (2013) Lazzati, D., Morsony, B. J., Margutti, R., & Begelman, M. C. 2013, The Astrophysical Journal, 765, 103
  * López-Cámara et al. (2014) López-Cámara, D., Morsony, B. J., & Lazzati, D. 2014, MNRAS, 442, 2202, doi: 10.1093/mnras/stu1016
  * Mignone et al. (2007) Mignone, A., Bodo, G., Massaglia, S., et al. 2007, ApJS, 170, 228, doi: 10.1086/513316
  * Mignone et al. (2012) Mignone, A., Zanni, C., Tzeferacos, P., et al. 2012, ApJS, 198, 7, doi: 10.1088/0067-0049/198/1/7
  * Parsotan (2021) Parsotan, T. 2021, parsotat/ProcessMCRaT: Cyclo-Synchrotron Release, v1.0.1, Zenodo, Zenodo, doi: 10.5281/zenodo.4918108
  * Parsotan et al. (2021) Parsotan, T., Cochrane, R. K., Hayward, C. C., et al. 2021, MNRAS, 501, 1591, doi: 10.1093/mnras/staa3765
  * Parsotan & Lazzati (2018) Parsotan, T., & Lazzati, D. 2018, ApJ, 853, 8, doi: 10.3847/1538-4357/aaa087
  * Parsotan & Lazzati (2021) Parsotan, T., & Lazzati, D. 2021, The Astrophysical Journal, 922, 257
  * Parsotan & Lazzati (2022) Parsotan, T., & Lazzati, D. 2022, ApJ, 926, 104, doi: 10.3847/1538-4357/ac4093
  * Parsotan et al. (2018) Parsotan, T., López-Cámara, D., & Lazzati, D. 2018, ApJ, 869, 103, doi: 10.3847/1538-4357/aaeed1
  * Pe’er (2008) Pe’er, A. 2008, ApJ, 682, 463, doi: 10.1086/588136
  * Pe’er & Ryde (2011) Pe’er, A., & Ryde, F. 2011, ApJ, 732, 49, doi: 10.1088/0004-637X/732/1/49
  * Pe’er & Ryde (2011) Pe’er, A., & Ryde, F. 2011, The Astrophysical Journal, 732, 49
  * Python Core Team (2019) Python Core Team. 2019, Python: A dynamic, open source programming language, Python Software Foundation.  https://www.python.org/
  * Rees & Meszaros (1994) Rees, M. J., & Meszaros, P. 1994, ApJ, 430, L93, doi: 10.1086/187446
  * Rees & Mészáros (2005) Rees, M. J., & Mészáros, P. 2005, ApJ, 628, 847, doi: 10.1086/430818
  * Stamatellos & Whitworth (2005) Stamatellos, D., & Whitworth, A. P. 2005, Astronomy & Astrophysics, 439, 153
  * Vurm & Beloborodov (2016) Vurm, I., & Beloborodov, A. M. 2016, ApJ, 831, 175, doi: 10.3847/0004-637X/831/2/175
  * Woosley & Heger (2006) Woosley, S., & Heger, A. 2006, The Astrophysical Journal, 637, 914
  * Yonetoku et al. (2004) Yonetoku, D., Murakami, T., Nakamura, T., et al. 2004, ApJ, 609, 935, doi: 10.1086/421285
  * Zhang & Yan (2011) Zhang, B., & Yan, H. 2011, ApJ, 726, 90, doi: 10.1088/0004-637X/726/2/90


# Network Analysis of the iNaturalist
Citizen Science Community

Yu Lu Liu and Thomas Jiralerspong  McGill UniversityMontrealCanada
<EMAIL_ADDRESS><EMAIL_ADDRESS>
(2022)

###### Abstract.

In recent years, citizen science has become a larger and larger part of the
scientific community. Its ability to crowd source data and expertise from
thousands of citizen scientists makes it invaluable. Despite the field’s
growing popularity, the interactions and structure of citizen science projects
are still poorly understood and under analyzed. We use the iNaturalist citizen
science platform as a case study to analyze the structure of citizen science
projects. We frame the data from iNaturalist as a bipartite network and use
visualizations as well as established network science techniques to gain
insights into the structure and interactions between users in citizen science
projects. Finally, we propose a novel unique benchmark for network science
research by using the iNaturalist data to create a network which has an
unusual structure relative to other common benchmark networks. We demonstrate
using a link prediction task that this network can be used to gain novel
insights into a variety of network science methods.

††copyright: acmcopyright††journalyear: 2022††doi:
XXXXXXX.XXXXXXX††conference: Workshop on Theories and Applications of Network
Science; December, 2022; Montreal, QC††price: 15.00††isbn:
978-1-4503-XXXX-X/18/06

Note to reviewers: since our project contributes a new graph dataset, we
structured the paper to be as intuitive as possible, so it does not follow the
proposed structure.

##  1\. Introduction and Motivation

The field of network science and the larger scientific community have long
benefited from one another. Network science methods have been used to analyze
relationships in the scientific community, such as co-citation, co-authorship,
etc. On the other side, datasets collected from the scientific community, for
instance the citation networks Cora, Citeseer and Pubmed (Sen et al., 2008),
are commonly used to train and to benchmark new network science methods. Most
works, however, are centered around academia: using and analysing connections
between academic papers or their authors. Communities of citizen scientists
are thus under-explored, both as a source of rich data and as an area of
application for network science methods.

Understanding citizen science communities is also an important task in its own
right. These communities contribute to science in two significant ways: i)
they democratize science by involving the public in the scientific process and
ii) they can be an important source of data in topics such as biodiversity,
astronomy, health, etc. As an example, Bombus irisanensis Cockerell, a
vulnerable bumblebee species which has not been reported to be observed since
the 1990’s, was ”rediscovered” in 2020 thanks to two observations published on
iNaturalist, a citizen science platform focused on gathering biodiversity
organisms (Wilson et al., 2020). Since citizen science platforms rely on the
participation of volunteers, they face challenges related to data quality
(cit, 2015), which is affected by the volunteers’ motivation, skill and
knowledge. The longevity of a platform also depends on recruitment, retention
and involvement of its participants. Understanding the communities that form
on these platforms are thus crucial in improving citizen science platforms,
both as a data collection tool and as a social network for citizen scientists.

Our main contributions are therefore:

  * •

We propose a novel collaboration network benchmark (built from a citizen
science community) to the network science community. This network has an
unusual structure compared to existing network science benchmarks.

  * •

We demonstrate using a link prediction task that this network’s unusual
structure allows it to be used to gain novel insights into network science
methods

  * •

We analyze this network using network science methods and draw insights on the
current state of the citizen science community we gathered data from.

##  2\. Related Work

There have been previous papers analyzing citizen science, but they do not
involve network science. For example, Aceves-Bueno et al. (Aceves-Bueno et
al., 2017) used statistical methods to analyze the accuracy of citizen science
data used in various publications. Bonnet et al. (Bonney et al., 2016) used
statistical methods to quantify the impact that citizen science has had on
scientific research, on the general public’s understanding of various
scientific domains, and on the personal well-being of people contributing to
citizen science projects. Other studies have specifically studied how citizen
scientists contribute to the iNaturalist platform using other techniques than
network science: Aristeidou et al. (Aristeidou et al., 2021) used various
statistical approaches to analyze the behavior of young people contributing to
the iNaturalist platform. Di Cecco et al. (Di Cecco et al., 2021) used various
statistical and clustering techniques to understand biases in citizen science
data and to understand which types of users use the iNaturalist platform. Our
work obviously differs from these because we propose to analyze citizen
science communities under the perspective of network science.

Finally, although combining network science and citizen science is a very
under-explored area, there have been a few previous studies in the area.
Pelacho et al. (Pelacho et al., 2020), for instance, used network science
techniques to analyze the evolution and collaboration networks of scientific
publications stemming from citizen science projects. Although this analysis is
related and applied to citizen science, it is still in the context of academic
publications and thus does not capture the full breadth of possibilities
provided by citizen scientists. Mazumdar et al. (Mazumdar and Thakker, 2020)
is an another example: they applied network science techniques in order to
analyze how people discuss citizen science on social media platforms such as
Twitter. This application of network science, while related to citizen
science, is not directly applied to citizen science communities but instead to
the way people communicate about citizen science.

To the best of our knowledge, our work is the first to use network science
techniques to analyze and understand the interactions within citizen science
communities in a non-academic setting, to view data from the iNaturalist
platform as a network, and to apply network science techniques to the
iNaturalist platform.

##  3\. Problem Definition

Our project focuses on citizen science communities on iNaturalist, a social
network platform where users share and map observations of nature. It is a
joint initiative of the California Academy of Sciences and the National
Geographic Society, launched in 2008. As of 2022, the number of observations
on the platform has reached 115 million (Wikipedia, 2022n). Users on the
platform can upload pictures of organisms they observe and identify organisms
in pictures that other users upload. Each observation can only reach ”research
grade” quality when there is a community consensus on the species of the
observed organism.

The first goal of our project is to gain insights about the communities on the
iNaturalist citizen science community:

  * •

Are there clusters of users who tend identify the same groups of species (e.g.
birdwatchers)?

  * •

What are the species that most often identified? Does their taxonomic class
come into play?

  * •

How do new users (e.g. users who created their accounts in 2022) fit into the
network?

  * •

etc.

To do this, in Section 5, we visualize and analyze the iNaturalist network to
identify interesting patterns that can be found in it.

The second goal of our project is to contribute a new benchmark network to
collaboration network research. To do this, in Section 6 we analyze
iNaturalist by viewing each taxon as a ”project” that users can contribute to
by correctly identifying observations of said taxon. Users that identify a
common taxon can thus be viewed as ”collaborators”, making the resulting
network a collaboration network.

In Section 6.2 we first investigate whether our ”collaboration” network, which
is built from a citizen science platform, has a different structure than
existing collaboration networks, which stem from academia. To do so, we
compute network measures commonly used to analyze co-authorship networks and
offer a comparison and an interpretation of the results.

In Section 6.3, we conduct further experiments to see if our network can
provide valuable insights into existing network science methods by checking
whether common link prediction methods perform differently or unusually on our
network compared to existing collaborating networks.

##  4\. Data Collection

Figure 1. Histogram of users’ account creation years. We see this information
is unavailable for around 30% of users.

To build the iNaturalist collaboration network, we start by collecting
observations. We collected all research-grade observations made in the Alberta
15th Division, a region situated within the Rocky Mountains in Western Canada
and includes the well-known Jasper and Banff national parks. The observations
are collected via GBIF 111 GBIF.org (27 November 2022) GBIF Occurrence
Download https://doi.org/10.15468/dl.wt7szw . This corresponds to a collection
of 41K observations identified by 2998 users and covering 1816 species. For
each observation, we log the species of the observed organism (as the
corresponding taxon ID) and the ID of the user who correctly identified it.

We define two types of nodes: a user node for each distinct user ID, and a
species node for each distinct taxon ID. We draw a link between a user node
and a species node if the user in question has correctly identified the
species at least once within our collected observations. We denote this
initial bipartite graph by user-taxon. By projecting the user-taxon network
onto the user nodes, we obtain the user-user network which contains only user
nodes. Two users are connected if they have ”collaborated”, i.e. identified at
least one common species within our collected observations. By projecting the
user-taxon network onto the taxon nodes, we obtain the taxon-taxon network
which contains only taxon nodes. Two taxons are connected if they share at
least one common user-identifier within our collected observations.

Figure 2. In A1 and A2, we visualize the bipartite graph user-taxon. In A1,
user nodes’ size is proportional to their (out) degree. In A2, taxon nodes’
size is proportional to their (in) degree. In B, we visualize the taxon-taxon
projection. In C, we visualize the user-user projection. User nodes’ colour
indicates whether they are new (joined on or after 2020). Taxon nodes’ colour
indicates their taxonomic class.

We add the following additional node attributes:

  * •

Taxonomic class (e.g. insects, birds, etc.) of each taxon to taxon nodes. This
additional information allows us to view similar species as a group.

  * •

Account creation year of each user to user nodes. This additional information
allows us to track new users and compare their involvement within the
community to older user accounts.

We note one limitation of our approach: the account creation year is an
unavailable information for some user accounts. This limits our capability to
differentiate between new users and old users, so we instead analyze users
known to be new (referred to as ”new users” from now on) versus other users.
We present a histogram of account creation years in Figure 1. We define new
users to be users with accounts created on or after 2020, which corresponds to
around 50% of the user base with a known account creation year.

##  5\. Analysis of the iNaturalist Community

Figure 3. The top 10 taxons with the highest number of user-identifiers (in-
degree). **We use ”solo identifiers” to refer to users who only identified the
taxon in question, within our collected observations. For instance, there are
80 users who identified bighorn sheeps and no other taxons. Source of
pictures: Wikipedia (Wikipedia, 2022e, i, j, l, d, c, g, h, r)

We present the visualizations of three networks in Figure 2: the bipartite
graph, the taxon-taxon projection and the user-user projection. The
visualization is produced using the Force Atlas algorithm via Gephi (Bastian
et al., 2009). Each user node is colour-coded to indicate whether the user is
new (joined after 2020) or not. Each taxon node is colour-coded to indicate
the taxonomic class to which it belongs:

  * •

Magnoliopsida, a class of flowering plants (Wikipedia, 2022p): 575 taxons;

  * •

Insecta, insects: 489 taxons;

  * •

Aves, birds (Wikipedia, 2022f): 215 taxons;

  * •

Agaricomycetes, a class of fungi which includes mushroom-forming fungi
(Wikipedia, 2022a): 180 taxons;

  * •

Liliopsida, lilies(Wikipedia, 2022o): 125 taxons;

  * •

Mammalia, mammals: 53 taxons;

  * •

Other, taxonomic classes not included above: 179 taxons.

Our main goal in this section is to draw insights about the community using
our networks.

###  5.1. user-taxon bipartite graph

In the bipartite graph (Figure 2: A1 and A2), we set edges to be directed from
users to taxons so that in-degree and out-degree values can better distinguish
between important user nodes and important taxon nodes.

First, we notice an obvious but small cluster at the bottom of the graph
consisting of an user who identified several Agaricomycetes. This user and
this taxonomic class seem to have little ”interaction” with the rest of the
network. If we exclude this small cluster, the bipartite graph has a peanut-
like shape: the left is largely occupied by birds, mammals and users who
identify species of those two classes; the right, on the other hand, is
largely occupied by taxons of the Magnoliopsida class and users who identify
these flowers. This suggests that the INaturalist community in the Alberta
15th Division can be viewed as one ”specialist” individual and two groups: one
focused on birds and mammals, and the other focused on Magnoliopsida flowers.

We can gain further insights into the community by visualizing the degree
values. In the first view of the graph (Figure 2: A1), the user node size is
proportional to the node’s out-degree value. In other words, big user nodes in
the visualization correspond to users who identified a larger amount of
distinct taxons. We can see that the biggest nodes are not new users and are
placed either within the right cluster or in the middle of the two main
clusters. This suggests that users who can identify a larger amount of species
tend not to be new users and that their focus is oriented more towards
Magnoliopsida flowers. Intuitively, this makes sense as Magnoliopsida is the
class containing the most taxons.

In the second view of the graph (Figure 2: A2), the taxon node size is
proportional to the node’s in-degree value. This means that a taxon node is
big if it has a large number of users who identified it. We can see that
although there are few mammal nodes, they are the biggest: there are more
users who identified mammals than taxons belonging to other taxonomic classes.
Moreover, we can observe that the left cluster tends to contain larger nodes
than the right cluster. This suggests that while there are fewer species of
mammals and birds, they are more easily identified by members of the
community. On the other hand, while there are many species of Magnoliopsida
flowers, they are in general harder to identify and are therefore able to be
identified by fewer users.

To discern the exact species that are identified by most users, we present in
Figure 3 the top 10 taxons with the highest in-degree value. We believe that
taxons on this list are either easily identifiable/observable (more
observations allow for more identifications), more interesting to
identify/observe, or very likely a combination of all these factors. We note
in particular that four of these taxons, namely the bighorn sheep (Wikipedia,
2022e), the golden-mantled ground squirrel (Wikipedia, 2022l), the Columbian
ground squirrel (Wikipedia, 2022j), and the Clark’s nutcracker (Wikipedia,
2022i), have the Western North American mountains as their one and only
habitat, which could explain their importance within the iNaturalist community
of the Alberta 15th Division.

Figure 4. Agaricomycetes taxons of interest. We note that most of these taxons
have distinctive colours or shape, which might make them easier or more
interesting to observe/identify. Source of pictures: Wikipedia
(wiki:Hydnellum_pecki; Wikipedia, 2022s, b, k, m, q)

We also highlight a group of users we call ”solo identifiers”, which are users
who only identified one taxon and no other taxons. If we sum up the number of
solo identifiers for the top 10 taxons presented in Figure 3, we obtain a
total of 326 such users, which corresponds to around 10% of the whole user
base. We elaborate on this group of users in Section 5.3.

###  5.2. taxon-taxon projection

In the taxon-taxon projection (Figure 2: B), the node size is proportional to
the node’s degree, which means that bigger nodes represent taxons that share
user-identifiers with a greater number of other taxons. Through this lense, we
can see that the majority of bigger nodes are clustered at the center,
surrounded by small clusters made of smaller nodes. By using the colour
coding, we can also observe that taxons of the same taxonomic class tend to be
grouped together, the most obvious cases being Agaricomycetes, Magnoliopsida,
and Aves (birds). This means that users tend to identify taxons belonging to
the same taxonomic class. For example, a user who identified a bird taxon
might also identify other bird taxons, which makes all bird taxons be grouped
closely together. Intuitively, this makes sense as nature enthusiasts tend to
be grouped as bird-watchers, mycophiles, etc.

Interestingly, we see bigger Agaricomycetes nodes that lie either within the
main cluster or between the main cluster and the far-away Agaricomycetes
cluster. We present their scientific name and their physical appearance in
Figure 4. Compared to Agaricomycetes nodes in the far-away cluster which are
identified by one ”specialist” individual, they might be more easily
identifiable, observable or more interesting to people, such that many users
who usually identify taxons of other taxonomic classes can identify them. This
small case study shows that the taxon-taxon projection can be useful for
gaining insights about the expertise of a citizen science community.

Figure 5. User-user projection with solo user-identifiers coloured in gray,
solo user-identifiers who identified one of the top 10 taxons (from 5.2) in
bright colors. Not that only 8 visible as there are 2 taxons containing too
few solo user-identifiers to be visible. Other users are coloured in black.

###  5.3. user-user projection

In the user-user projection (Figure 2: C), we also set the node size to be
proportional to the node’s degree. In this case, bigger nodes represent users
who share a greater number of ”projects” with other users, each ”project”
being the identification of a taxon. We see that bigger nodes (more
collaborative users) are mostly in the middle while smaller nodes (less
collaborative users) are pushed to the periphery and form small clusters. To
further investigate these clusters, we return to the ”solo identifiers”
introduced in Section 5.1. In Figure 5, we colour-code the solo identifiers so
that their user nodes are either black, gray, or brightly coloured. Each
bright colour denotes solo-identifiers of one of the top 10 taxons listed in
Section 5.1, whereas gray denotes solo-identifiers of other taxons and black
denotes other users. We see that each small cluster corresponds exactly to the
solo-identifiers of one taxon and that all solo-identifier nodes are pushed to
the periphery. This shows that the visualization of the user-user projection
can be used to identify groups of users who behave in a significantly
different way than others.

In Figure 2 (C) we also colour-code the new users, but we do not see any
obvious trend.

|  iNaturalist (ours) |  GR-QC |  Erdos |  NetScience |  CS-PhD
---|---|---|---|---|---
Number of nodes |  2998 |  4158 |  5094 |  379 |  1882
Number of edges |  585115 |  13422 |  7515 |  914 |  1740
Number of connected components |  3 |  1 |  17 |  1 |  168
Size of largest component |  2994 |  4158 |  4991 |  379 |  1025
Power law fit slope |  -0.475 |  -1.974 |  -1.678 |  -1.553 |  -2.136
Average shortest path length |  1.906 |  6.049 |  5.511 |  6.042 |  11.659
Average clustering coefficient |  0.806 |  0.557 |  0.0818 |  0.741 |  0.00507
Spectral gap |  3818.504 |  34.999 |  28.143 |  6.431 |  5.969
Degree vs. degree correlation coefficient |  -0.27 |  0.639 |  -0.444 |  -0.0817 |  -0.197
Degree vs. clustering coefficient correlation coefficient |  -0.745 |  0.0513 |  0.0153 |  -0.362 |  0.0195
Table 1. Basic metrics for all collaboration networks. Outliers for the
iNaturalist network are in bold

##  6\. Comparison to other collaboration networks

As mentioned previously, the user-user network can be viewed as a
collaboration network where users are connected if they contributed to the
identification ”project” of a common species. We thus compare this network to
existing collaboration graphs from academic scientific communities. We conduct
this comparison first by computing common network metrics and comparing their
values for our network and for the academic collaboration networks to see if
our network has a distinct structure compared to thee academic networks. We
then compare the performance of several common link prediction methods on our
network compared to the academic collaboration networks to see if our network
can provide interesting insights into the performance of these methods.

###  6.1. Other Collaboration Networks

We compare our user-user projection network to 4 other benchmark collaboration
networks of similar size which we found in the Network Repository (Rossi and
Ahmed, 2015):

#### GR-QC

(Leskovec et al., 2007)

The most common form of collaboration in the academic world is co-authorship.
The GR-QC collaboration network contains information about scientists who have
published papers in the fields of general relativity and quantum cosmology.
Each node is a scientist and there is an edge between two nodes $u$ and $v$ if
scientist $u$ has co-authored a paper with scientist $v$.

#### Erdos

(Batagelj and Mrvar, 2000)

Paul Erdos was one of the most productive mathematicians in history,
publishing more than 1500 papers during his lifetime. The Erdos collaboration
network is a network of Erdos, his collaborators, and his collaborators’
collaborators. The network is defined in the same way as GR-QC.

#### NetScience

(Newman, 2006)

The NetScience collaboration network contains information about 379 scientists
who have conducted research in the field of network science. The network is
defined in the same way as GR-QC.

#### CS-PhD

(De Nooy et al., 2011)

Although collaboration networks are often studied in the context of co-
authorship, another collaborative relationship which exists in the academic
world is the PhD-advisor relationship.

In the CS-PhD collaboration network, each node is a person who holds a PhD in
theoretical computer science and there is an edge between two nodes $u$ and
$v$ if person $u$ was the PhD advisor for person $v$.

###  6.2. Basic Metrics

To compare our collaboration network to other academic collaboration networks,
we first compute some basic metrics for our network as well as the benchmark
networks described in section 6.1. The results are displayed in Table 1, with
outliers for the iNaturalist network in bold. We now analyze these outliers.

The iNaturalist network has a much larger number of edges compared to the
other networks, despite having a similar number of nodes. This shows that the
iNaturalist network is much denser than the academic networks, indicating that
collaborations are much more frequent in citizen science compared to academic
science. This is expected because participating in a citizen science taxon
”project” is much more accessible than collaborating academically.

The iNaturalist network also has a much smaller average shortest path length
compared to the academic networks. This shows that the iNaturalist network is
much more small-world compared to the academic networks, which is expected due
to the iNaturalist network having many more edges than the academic networks.

Finally, the iNaturalist network has a much larger spectral gap compared to
the other networks. A higher spectral gap indicates that a graph is more
connected, so this quantitatively shows that the iNaturalist network is much
more connected than the academic science networks.

These three metrics all indicate that citizen scientists tend to collaborate
much more often than academic scientists. This confirms the intuitive notion
that, although academic scientists may be more qualified than citizen
scientists, citizen science compensates by making it easier to collaborate
with others, thus allowing for a much larger number of collaborations.

We can also see from these metrics that the iNaturalist network has a very
different structure compared to existing academic collaboration network
benchmarks, indicating that network science algorithms may exhibit very
different behavior when evaluated on the iNaturalist network compared to
academic collaboration network. This might make it a valuable benchmark for
network scientists to test if their algorithms are robust to many different
types of network structures. We test this hypothesis using a link prediction
task in the following section.

###  6.3. Link Prediction in user-user

|  iNaturalist (ours) |  GR-QC |  Erdos |  NetScience |  CS-PhD
---|---|---|---|---|---
Jaccard Coefficient |  0.918 |  0.903 |  0.620 |  0.906 |  0.503
Resource Allocation Index |  0.985 |  0.903 |  0.620 |  0.910 |  0.503
Adamic Adar Index |  0.960 |  0.903 |  0.620 |  0.910 |  0.503
Preferential Attachment |  0.915 |  0.739 |  0.545 |  0.638 |  0.345
DeepWalk |  0.843 |  0.681 |  0.857 |  0.640 |  0.647
Table 2. ROC-AUC scores for different link prediction methods on different
collaboration network datasets

We evaluate 5 link prediction methods (Jaccard Coefficient, Resource
Allocation Index, Adamic Adar Index, Preferential Attachment, DeepWalk
(Perozzi et al., 2014)) on our collaboration network as well as the 4
collaboration networks from section 6.1. To do this we first randomly delete
20% of the edges in the network. We call these edges ”gold edges” and the
edges which remain in the network ”remaining edges”. We then create a list of
”candidate edges” by generating a complete graph on all the nodes in the
network and deleting from it the edges which are part of the ”remaining
edges”.

Finally we frame link prediction as a classification task where our link
prediction method must predict whether each of the candidate edges is a gold
edge or not while only having access to the remaining edges. We compute the
ROC-AUC scores for each method and present the results in Table 2.

We see that almost all the methods, with the exception of DeepWalk, perform
the best on the iNaturalist network compared to the other networks, and that
even DeepWalk performs second best on it. This suggests that it is easier to
predict which citizen scientists have contributed to at least one of the same
”projects” than it is to predict which academic scientists have co-authored a
paper together.

We also see that, for the academic networks, the 3 first methods (Jaccard
Coefficient, Resource Allocation Index, Adamic Adar Index) all perform almost
identically. It is thus very hard to compare these methods or say which one of
them is better using these academic networks. With the iNaturalist user-user
network, however, we see that the 3 methods’ performances differ more from
each other, thus allowing us to compare them and rank them.

We also see that the 4th method (Preferential Attachment) performs much worse
than the first 3 methods on the 4 academic networks, whereas it has almost
identical performance to Jaccard Coefficient and only slightly worst
performance than Resource Allocation Index and Adamic Adar Index on
iNaturalist. This suggests that Preferential Attachment performs better on
denser graphs.

DeepWalk seems to be an outlier method as its performance on the various
networks does not follow the trend of the other methods: it has the worst
performance on GR-QC, NetScience, and our dataset, but has the best
performance on Erdos and CS-PhD. We hypothesize that this might be due to the
fact that DeepWalk is the only GNN based method that we tested. GNN based
methods are known to be very sensitive to the hyperparameters used to train
them and we regrettably did not have time to do a proper hyperparameter search
for DeepWalk. This may be why its results are outliers compared to the other
methods.

These insights derived from the methods’ performance on our network show that
the iNaturalist use-user network is a valuable valuable benchmark for
differentiating between methods when other benchmarks are unable to do so, and
for gaining new insights into how certain methods achieve better performance
on certain types of networks.

##  7\. Conclusion

In this paper we used network science to gain insights from the iNaturalist
community in the Alberta 15th Division. Through various visualizations and
degree measures, we found that:

  * •

The INaturalist community in the Alberta 15th Division can be viewed as one
individual specializing in Agaricomycetes, and two groups: one focused on
birds and mammals, and the other focused on Magnoliopsida flowers.

  * •

Users who can identify a larger amount of species tend not to be new users and
that their focus is oriented more towards Magnoliopsida flowers.

  * •

Users tend to identify taxons belonging to the same taxonomic class.

Most importantly, we demonstrated how the bipartite network and its two
projection networks can be used to analyze a citizen science community.

We also introduced the iNaturalist collaboration network as a new benchmark
for network science research. We showed that it has an unusual structure
compared to existing collaboration network benchmarks and demonstrated using a
link prediction task that this unusual structure makes it a valuable benchmark
for evaluating and analyzing the properties of network science methods.

##  8\. Future Work

For future work, we would like to obtain more data from the iNaturalist
platform to be able to build even larger benchmark data sets. We would also
like to analyze how different community detection methods perform on detecting
the different taxonomic classes in the taxon-taxon projection network.

##  9\. Acknowledgements

We would like to thank Professor Reihaneh Rabbany as well as all our reviewers
for their valuable suggestions and feedback all throughout this project.

## References

  * (1)
  * cit (2015) 2015\.  Rise of the citizen scientist.  _Nature_ 524, 7565 (01 Aug 2015), 265–265.  https://doi.org/10.1038/524265a
  * Aceves-Bueno et al. (2017) Eréndira Aceves-Bueno, Adeyemi S Adeleye, Marina Feraud, Yuxiong Huang, Mengya Tao, Yi Yang, and Sarah E Anderson. 2017\.  The accuracy of citizen science data: a quantitative review.  _Bulletin of the Ecological Society of America_ 98, 4 (2017), 278–290.
  * Aristeidou et al. (2021) Maria Aristeidou, Christothea Herodotou, Heidi L Ballard, Alison N Young, Annie E Miller, Lila Higgins, and Rebecca F Johnson. 2021\.  Exploring the participation of young citizen scientists in scientific research: The case of iNaturalist.  _Plos one_ 16, 1 (2021), e0245682.
  * Bastian et al. (2009) Mathieu Bastian, Sebastien Heymann, and Mathieu Jacomy. 2009\.  Gephi: An Open Source Software for Exploring and Manipulating Networks.  (2009).  http://www.aaai.org/ocs/index.php/ICWSM/09/paper/view/154
  * Batagelj and Mrvar (2000) Vladimir Batagelj and Andrej Mrvar. 2000.  Some analyses of Erdos collaboration graph.  _Social Networks_ 22, 2 (2000), 173–186.
  * Bonney et al. (2016) Rick Bonney, Tina B Phillips, Heidi L Ballard, and Jody W Enck. 2016.  Can citizen science enhance public understanding of science?  _Public understanding of science_ 25, 1 (2016), 2–16.
  * De Nooy et al. (2011) Wouter De Nooy, Andrej Mrvar, and Vladimir Batagelj. 2011\.  _Exploratory social network analysis with Pajek_. Vol. 27.  Cambridge University Press.
  * Di Cecco et al. (2021) Grace J Di Cecco, Vijay Barve, Michael W Belitz, Brian J Stucky, Robert P Guralnick, and Allen H Hurlbert. 2021.  Observing the observers: How participants contribute data to iNaturalist and implications for biodiversity science.  _BioScience_ 71, 11 (2021), 1179–1188.
  * Leskovec et al. (2007) Jure Leskovec, Jon Kleinberg, and Christos Faloutsos. 2007\.  Graph Evolution: Densification and Shrinking Diameters.  _ACM Trans. Knowl. Discov. Data_ 1, 1 (mar 2007), 2–es.  https://doi.org/10.1145/1217299.1217301
  * Mazumdar and Thakker (2020) Suvodeep Mazumdar and Dhavalkumar Thakker. 2020.  Citizen Science on Twitter: Using Data Analytics to Understand Conversations and Networks.  _Future Internet_ 12, 12 (2020).  https://doi.org/10.3390/fi12120210
  * Newman (2006) Mark EJ Newman. 2006\.  Finding community structure in networks using the eigenvectors of matrices.  _Physical review E_ 74, 3 (2006), 036104.
  * Pelacho et al. (2020) M. Pelacho, G. Ruiz, F. Sanz, A. Tarancon, and J. Clemente-Gallardo. 2020.  Analysis of the evolution and collaboration networks of citizen scientific publications.  _AI Magazine_ 126, 3 (Oct. 2020), 225\.  https://doi.org/10.1007/s11192-020-03724-x
  * Perozzi et al. (2014) Bryan Perozzi, Rami Al-Rfou, and Steven Skiena. 2014\.  DeepWalk. In _Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining_. ACM.  https://doi.org/10.1145/2623330.2623732
  * Rossi and Ahmed (2015) Ryan A. Rossi and Nesreen K. Ahmed. 2015.  The Network Data Repository with Interactive Graph Analytics and Visualization. In _AAAI_.  https://networkrepository.com
  * Sen et al. (2008) Prithviraj Sen, Galileo Namata, Mustafa Bilgic, Lise Getoor, Brian Galligher, and Tina Eliassi-Rad. 2008\.  Collective Classification in Network Data.  _AI Magazine_ 29, 3 (Sep. 2008), 93\.  https://doi.org/10.1609/aimag.v29i3.2157
  * Wikipedia (2022a) Wikipedia. 2022a.  Agaricomycetes — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Agaricomycetes&oldid=1113445436.  [Online; accessed 05-December-2022].
  * Wikipedia (2022b) Wikipedia. 2022b.  Alloclavaria purpurea — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Alloclavaria%20purpurea&oldid=1051410601.  [Online; accessed 05-December-2022].
  * Wikipedia (2022c) Wikipedia. 2022c.  American black bear — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=American%20black%20bear&oldid=1121916752.  [Online; accessed 05-December-2022].
  * Wikipedia (2022d) Wikipedia. 2022d.  Aquilegia flavescens — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Aquilegia%20flavescens&oldid=1094178920.  [Online; accessed 05-December-2022].
  * Wikipedia (2022e) Wikipedia. 2022e.  Bighorn sheep — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Bighorn%20sheep&oldid=1125781182.  [Online; accessed 05-December-2022].
  * Wikipedia (2022f) Wikipedia. 2022f.  Bird — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Bird&oldid=1124903835.  [Online; accessed 05-December-2022].
  * Wikipedia (2022g) Wikipedia. 2022g.  Canada jay — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Canada%20jay&oldid=1122059444.  [Online; accessed 05-December-2022].
  * Wikipedia (2022h) Wikipedia. 2022h.  Chamaenerion angustifolium — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Chamaenerion%20angustifolium&oldid=1125480878.  [Online; accessed 05-December-2022].
  * Wikipedia (2022i) Wikipedia. 2022i.  Clark’s nutcracker — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Clark’s%20nutcracker&oldid=1122867502.  [Online; accessed 05-December-2022].
  * Wikipedia (2022j) Wikipedia. 2022j.  Columbian ground squirrel — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Columbian%20ground%20squirrel&oldid=1124139650.  [Online; accessed 05-December-2022].
  * Wikipedia (2022k) Wikipedia. 2022k.  Coprinus comatus — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Coprinus%20comatus&oldid=1106032411.  [Online; accessed 05-December-2022].
  * Wikipedia (2022l) Wikipedia. 2022l.  Golden-mantled ground squirrel — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Golden-mantled%20ground%20squirrel&oldid=1121777526.  [Online; accessed 05-December-2022].
  * Wikipedia (2022m) Wikipedia. 2022m.  Hydnellum caeruleum — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Hydnellum%20caeruleum&oldid=1110709335.  [Online; accessed 05-December-2022].
  * Wikipedia (2022n) Wikipedia. 2022n.  INaturalist — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=INaturalist&oldid=1122751314.  [Online; accessed 05-December-2022].
  * Wikipedia (2022o) Wikipedia. 2022o.  Liliopsida — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Liliopsida&oldid=1062713963.  [Online; accessed 05-December-2022].
  * Wikipedia (2022p) Wikipedia. 2022p.  Magnoliopsida — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Magnoliopsida&oldid=1125274045.  [Online; accessed 05-December-2022].
  * Wikipedia (2022q) Wikipedia. 2022q.  Phellinus tremulae — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Phellinus%20tremulae&oldid=1008806146.  [Online; accessed 05-December-2022].
  * Wikipedia (2022r) Wikipedia. 2022r.  Red deer — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Red%20deer&oldid=1125200239.  [Online; accessed 05-December-2022].
  * Wikipedia (2022s) Wikipedia. 2022s.  Sarcodon imbricatus — Wikipedia, The Free Encyclopedia.  http://en.wikipedia.org/w/index.php?title=Sarcodon%20imbricatus&oldid=1121756904.  [Online; accessed 05-December-2022].
  * Wilson et al. (2020) Joseph S. Wilson, Aaron D. Pan, David Emmanuel M. General, and Jonathan B. Koch. 2020.  More eyes on the prize: an observation of a very rare, threatened species of Philippine Bumble bee, Bombus irisanensis, on iNaturalist and the importance of citizen science in conservation biology.  _Journal of Insect Conservation_ 24, 4 (01 Aug 2020), 727–729.  https://doi.org/10.1007/s10841-020-00233-3