\documentclass[10pt]{article}

\usepackage{fullpage}
\usepackage{setspace}
\usepackage{parskip}
\usepackage{titlesec}
\usepackage[section]{placeins}
\usepackage{xcolor}
\usepackage{breakcites}
\usepackage{lineno}
\usepackage{hyphenat}





\PassOptionsToPackage{hyphens}{url}
\usepackage[colorlinks = true,
            linkcolor = blue,
            urlcolor  = blue,
            citecolor = blue,
            anchorcolor = blue]{hyperref}
\usepackage{etoolbox}
\makeatletter
\patchcmd\@combinedblfloats{\box\@outputbox}{\unvbox\@outputbox}{}{%
  \errmessage{\noexpand\@combinedblfloats could not be patched}%
}%
\makeatother


\usepackage{natbib}




\renewenvironment{abstract}
  {{\bfseries\noindent{\abstractname}\par\nobreak}\footnotesize}
  {\bigskip}

\titlespacing{\section}{0pt}{*3}{*1}
\titlespacing{\subsection}{0pt}{*2}{*0.5}
\titlespacing{\subsubsection}{0pt}{*1.5}{0pt}


\usepackage{authblk}


\usepackage{graphicx}
\usepackage[space]{grffile}
\usepackage{latexsym}
\usepackage{textcomp}
\usepackage{longtable}
\usepackage{tabulary}
\usepackage{booktabs,array,multirow}
\usepackage{amsfonts,amsmath,amssymb}
\providecommand\citet{\cite}
\providecommand\citep{\cite}
\providecommand\citealt{\cite}
% You can conditionalize code for latexml or normal latex using this.
\newif\iflatexml\latexmlfalse
\AtBeginDocument{\DeclareGraphicsExtensions{.pdf,.PDF,.eps,.EPS,.png,.PNG,.tif,.TIF,.jpg,.JPG,.jpeg,.JPEG}}





\usepackage{float}
















\iflatexml
\documentclass{article}
\fi

  % \usepackage[T1]{fontenc} % Not (yet) supported by LaTeXML.

  \usepackage[utf8]{inputenc}

% \usepackage[english]{babel} % Not (yet) supported by LaTeXML.

\usepackage{comment}

\usepackage[a4paper, total={6in, 8in}]{geometry}

\usepackage{amsmath}

% \usepackage{a4wide} % Not (yet) supported by LaTeXML.

\newcommand{\me}[][]{\mathrm{e}}

\usepackage{amsfonts}

% \usepackage{mathtools} % Not (yet) supported by LaTeXML.

\usepackage{enumerate}

\usepackage{rotating}

\usepackage{float}

\newcommand{\lv}[][]{\langle}

\newcommand{\rv}[][]{\rangle}

\usepackage{amsthm}

\usepackage{amssymb}

\usepackage{float}

\usepackage{natbib}

% \usepackage{authblk} % Not (yet) supported by LaTeXML.

\usepackage{tabularx}

\usepackage{multirow}

\usepackage{multicol}

\usepackage[table, dvipsnames]{xcolor}

% \usepackage{tablestyles} % Not (yet) supported by LaTeXML.

\usepackage{longtable}
% \usepackage{sectsty} % Not (yet) supported by LaTeXML.

\usepackage{graphicx}

\usepackage{amsmath}

\usepackage{amssymb}

\usepackage{multicol}

% \usepackage{changes} % Not (yet) supported by LaTeXML.

\usepackage{bbm}

\usepackage{float}

\usepackage{booktabs}

% \usepackage{subcaption} % Not (yet) supported by LaTeXML.

\usepackage{rotating}

\usepackage[font={small,it}]{caption}

\usepackage{multirow}

\usepackage{tabularx}

\usepackage{natbib}

% \usepackage{authblk} % Not (yet) supported by LaTeXML.

\usepackage[colorlinks=true, linkcolor=black, citecolor=black, urlcolor=black]{hyperref}
% \usepackage{svg} % Not (yet) supported by LaTeXML.
\newcommand*\Laplace[][]{\mathop{}\!\mathbin\bigtriangleup}

\newcommand*\DAlambert[][]{\mathop{}\!\mathbin\Box}

\setlength\parindent{20pt}
\usepackage{setspace}

\usepackage{lineno}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{corollary}{Corollary}[theorem]
\newtheorem{lemma}[theorem]{Lemma}
\theoremstyle{definition}
\newtheorem{definition}{Definition}[section]
\theoremstyle{remark}
% First macro on next line not (yet) supported by LaTeXML: 
% \newtheorem*{remark}{Remark}
\newcommand{\overbar}[1][]{\mkern 1.5mu\overline{\mkern-1.5mu#1\mkern-1.5mu}\mkern 1.5mu}


\begin{document}

\title{A structural model for the process of collecting biodiversity data}



\author[1,2]{Kwaku P. Adjei}%
\author[3]{Jorge Sicacha-Parada}%
\author[1]{Ingelin Steinsland}%
\author[1,2]{Robert B. O'Hara}%
\affil[1]{Department of Mathematical Sciences, Norwegian University of Science and Technology, Trondheim - Norway}%
\affil[2]{Center for Biodiversity Dynamics, Norwegian University of Science and Technology, Trondheim - Norway.}%
\affil[3]{Norwegian Computing Center, Oslo - Norway}%


\vspace{-1em}



  
  \date{June 14, 2024}


\begingroup
\let\center\flushleft
\let\endcenter\endflushleft
\maketitle
\endgroup








\sloppy


\section*{Introduction}
Biodiversity data from surveys and other monitoring programs such as Citizen Science (called `CS' hereafter) play a vital role in estimating species distributions and making conservation decisions  \hyperref[csl:1]{(Bird et al., 2014)}. With an influx of data available to researchers through the various biodiversity databases such as iNaturalist \hyperref[csl:2]{(Matheson, 2014)}, GBIF \hyperref[csl:3]{(Telenius, 2011)}, Artsdatabanken (\url{https://www.artsdatabanken.no/}), eBird \hyperref[csl:4]{(Sullivan et al., 2014)}, amongst many others, there has been increased coverage and inferences about species distributions. However, these biodiversity data (especially those from CS) are subject to sampling and systematic biases, raising concerns about their use in scientific research \hyperref[csl:5]{(Lewandowski \& Specht, 2015}; \hyperref[csl:6]{Crall et al., 2011)}. The sampling bias may be due to sampling variation in space and time, and the systematic bias could arise because of species misreporting and misidentification as well as imperfect detection, amongst many others \hyperref[csl:1]{(Bird et al., 2014}; \hyperref[csl:7]{Bal{\'a}zs et al., 2021)}.

As discussed by \hyperref[csl:1]{(Bird et al., 2014)}, various statistical approaches have been used in analysing biodiversity data such as the generalised linear (mixed) models \hyperref[csl:6]{(Crall et al., 2011}; \hyperref[csl:8]{Delaney et al., 2008}; \hyperref[csl:9]{Brunsdon \& Comber, 2012)}, hierarchical models such as N-Mixture models (\hyperref[csl:10]{(K{\'e}ry \& Royle, 2015}; \hyperref[csl:11]{K{\'e}ry \& Royle, 2020)} and references therein) and joint species distribution models (\hyperref[csl:12]{(Ovaskainen \& Abrego, 2020)} and references therein). These methods depend on the type of response variable. Some modellers of biodiversity data use Maxent \hyperref[csl:13]{(Phillips et al., 2006}; \hyperref[csl:14]{Phillips et al., 2009)},  a constrained optimisation approach that finds the optimal species density subject to constraints \hyperref[csl:15]{(Chakraborty et al., 2011)}. Attaching any uncertainty to Maxent's predictions \hyperref[csl:16]{(Gelfand \& Shirota, 2019)} is impossible as a non-stochastic approach. However, \hyperref[csl:17]{(Renner et al., 2015)} has shown that Maxent is equivalent to point process models and that the point process framework can be used to obtain standard errors on model coefficients and predictions. Biodiversity data have also been modelled as typical geostatistic data with a binary response by including pseudo-absences \hyperref[csl:18]{(Ferrier et al., 2002}; \hyperref[csl:19]{Barbet-Massin et al., 2012)}. These pseudo-absences are used to account for the biases in Maxent to generate reliable and unbiased species distribution models 
 \hyperref[csl:20]{(Phillips et al., 2009}; \hyperref[csl:21]{Syfert et al., 2013}; \hyperref[csl:19]{Barbet-Massin et al., 2012}; \hyperref[csl:22]{Zaniewski et al., 2002)}. However, this approach adds an arbitrary amount of data and ignores the spatial autocorrelation between absences, as pointed out in \hyperref[csl:16]{(Gelfand \& Shirota, 2019)}. However, some studies have explored different methods to control the selection of pseudo-absences to fit the species distribution models  \hyperref[csl:20]{(Phillips et al., 2009}; \hyperref[csl:19]{Barbet-Massin et al., 2012}; \hyperref[csl:22]{Zaniewski et al., 2002}; \hyperref[csl:23]{L{\"u}tolf et al., 2006}; \hyperref[csl:24]{Engler et al., 2004}; \hyperref[csl:25]{Chefaoui \& Lobo, 2008)}.

Furthermore, other approaches propose modelling biodiversity data as a thinned point pattern \hyperref[csl:15]{(Chakraborty et al., 2011}; \hyperref[csl:26]{Fithian et al., 2015}; \hyperref[csl:27]{Sicacha-Parada et al., 2021)}. A point pattern is a collection of points with random locations. That is, the locations of the points are not fixed or previously chosen. A spatial point process is the model that determines the location and amount of points in an area. Thinning is an operation on point patterns that uses a specified rule to determine which points in the point pattern are deleted. For biodiversity data, the thinning of the point pattern of the true species occurrences is caused by inherent biases in them. For instance,  \hyperref[csl:27]{(Sicacha-Parada et al., 2021)} proposed that a single covariate related to a source of bias (e.g. accessibility) should be included in the linear predictor as a log-linear function to determine the probability of retaining an occurrence (i.e. a variable restricted to $(-\infty,0])$ as described in \hyperref[csl:28]{(Yuan et al., 2017)}). \hyperref[csl:26]{(Fithian et al., 2015)} and \hyperref[csl:29]{(Geldmann et al., 2016)}, however, proposed integrating extra sources of information, such as professional surveys, to account for these biases. These biases can be explicitly modelled as a function of known covariates if that information is available or as an extra random effect when data related to the biases is unavailable \hyperref[csl:30]{(Simmonds et al., 2020)}. However, these approaches do not consider how the various sources of biases discussed in the next paragraphs affect the observed data.

Firstly, biodiversity data are typically affected by biased sampling and unbiased imperfect sampling processes. Observers, mostly from CS projects, usually choose where they go. This decision is frequently influenced by factors such as accessibility \hyperref[csl:31]{(Monsarrat et al., 2019)} and where observers expect to find more occurrences (i.e. preferential sampling; \hyperref[csl:32]{(Diggle et al., 2010)}). This is known as the biased sampling process. On the other hand, some areas will be inaccessible due to government regulations or because they are protected areas \hyperref[csl:33]{(Inman et al., 2021)} or a sub-sample of biological materials are taken at a given location \hyperref[csl:34]{(Shimadzu et al., 2016)}. This sampling incompleteness is an example of an unbiased imperfect sampling process. \hyperref[csl:15]{(Chakraborty et al., 2011}; \hyperref[csl:30]{Simmonds et al., 2020)} and \hyperref[csl:27]{(Sicacha-Parada et al., 2021)} have acknowledged the role of the sampling process as a thinning factor in the context of point patterns. \hyperref[csl:26]{(Fithian et al., 2015)} explored a model that assumes the thinning process as a log-linear function that affects the intensity of the observed point pattern. Further work by \hyperref[csl:15]{(Chakraborty et al., 2011}; \hyperref[csl:27]{Sicacha-Parada et al., 2021}; \hyperref[csl:35]{Ahmad Suhaimi et al., 2021)} and \hyperref[csl:30]{(Simmonds et al., 2020)} explored the implications of not properly accounting for the sampling bias in biodiversity data. This was done by comparing improvements in goodness-of-fit and ecological interpretability of a model that accounts for sampling bias to approaches such as Maxent or a model that does not account for the sampling bias at all. In both cases, modelling the sampling bias as a thinning operation on a point pattern reduced the bias in the estimated effect of ecological covariates on the spatial distribution of species occurrences. However, these developments have not explored other sources of bias that can affect the process of collecting biodiversity data.

Secondly, imperfect detection is also another source of bias in biodiversity data \hyperref[csl:36]{(Mugumaarhahama et al., 2022}; \hyperref[csl:37]{K{\'e}ry \& Schmid, 2004}; \hyperref[csl:38]{Metz et al., 2020}; \hyperref[csl:10]{K{\'e}ry \& Royle, 2015}; \hyperref[csl:39]{MacKenzie et al., 2002)}. Imperfect detection occurs when the observer fails to detect the species even though it was present. Detectability must be accounted for in estimating species trend and abundance because some species are overlooked, and the detection of the species will depend on its behaviour \hyperref[csl:37]{(K{\'e}ry \& Schmid, 2004)}. The detection and identification of these species are influenced by the observer's attention to a particular species and place, time and factors determining visibility, such as weather conditions \hyperref[csl:40]{(Arazy \& Malkinson, 2022)}. Failure to account for the imperfect detection when analysing biodiversity data may bias parameter estimates \hyperref[csl:41]{(Welsh et al., 2013)} and thus make statistical inference less reliable \hyperref[csl:36]{(Mugumaarhahama et al., 2022}; \hyperref[csl:42]{Royle \& Link, 2006)}. Integrated species distribution models (hereafter, iSDMs) have been proposed to analyse these observed biodiversity data from various sources and sampling protocols \hyperref[csl:43]{(Koshkina et al., 2017}; \hyperref[csl:44]{Dorazio, 2014}; \hyperref[csl:45]{Erickson \& Smith, 2021)}. These iSDMs capitalize on the strengths of each dataset to better capture species distributions and dynamics \hyperref[csl:46]{(Isaac et al., 2020)}.

During most biodiversity data collection, more than one species are observed or reported. The report of one species can be misclassified as another. For example, iNaturalist reports potential misclassified species when a species of interest is queried on their website (\url{www.inaturalist.org}). Therefore, it is essential that we treat biodiversity data as multispecies data with possible misclassifications and model it as such \hyperref[csl:47]{(Adjei et al., 2023)}. Misidentification and/or misreporting species and other sources of false positives (collectively known as misclassification)  are critical issues to consider in biodiversity data.  Various methods have been developed to model these false positives, such as those proposed by \hyperref[csl:48]{(Miller et al., 2011}; \hyperref[csl:49]{Chambert et al., 2015}; \hyperref[csl:50]{Wright et al., 2020}; \hyperref[csl:51]{Strickfaden et al., 2020}; \hyperref[csl:11]{K{\'e}ry \& Royle, 2020}; \hyperref[csl:47]{Adjei et al., 2023)} and \hyperref[csl:10]{(K{\'e}ry \& Royle, 2015)} 
and references therein. These methods can be model-based, which includes taking a subset of the data with verifiable certainty, instructing observers only to record observations they are sure about, and increasing observer experiences \hyperref[csl:51]{(Strickfaden et al., 2020)}. These methods can also be design-based, for example, dependent double-observer method \hyperref[csl:51]{(Strickfaden et al., 2020)}. Failure to account for misclassification in biodiversity data modelling can increase bias and decrease the precision of the parameter estimates \hyperref[csl:52]{(Royle et al., 2007}; \hyperref[csl:1]{Bird et al., 2014}; \hyperref[csl:51]{Strickfaden et al., 2020)}, leading to accidental culling of endangered species and assessments of the population status and incorrect conservation decisions \hyperref[csl:53]{(Austen et al., 2016)}.

Accounting for the biases existing in biodiversity data is paramount for users of these data. However, there is no consensus on how this can be done, and it is indeed a growing research field within both statistics and ecology \hyperref[csl:54]{(Isaac et al., 2014}; \hyperref[csl:27]{Sicacha-Parada et al., 2021}; \hyperref[csl:47]{Adjei et al., 2023}; \hyperref[csl:55]{Magurran \& McGill, 2010}; \hyperref[csl:56]{Ahrends, 2010)}. In this paper, we propose a hierarchical multispecies model for biodiversity data that simultaneously accounts for multiple sources of bias and specifies each of them probabilistically. Some work has been done to propose a framework for observer-based biases in CS data \hyperref[csl:57]{(Arazy \& Malkinson, 2021)}, but no work has been done to develop a Bayesian framework that provides insight into modelling biases in biodiversity data. It relies on a straightforward specification of the observed biodiversity data as a thinned point pattern. This point pattern is affected by common biases such as uneven sampling effort, imperfect detection and misclassification. We propose three stages of thinning, one for each type of bias considered. Each stage produces a probability of retaining a point from the previous stages. This novel approach is flexible and can accommodate more biases in biodiversity data beyond the ones discussed in this paper. This paper aims to introduce our modelling framework and its properties.

The proposed framework in this paper is not excused from identifiability issues \hyperref[csl:44]{(Dorazio, 2014}; \hyperref[csl:43]{Koshkina et al., 2017}; \hyperref[csl:26]{Fithian et al., 2015)}. A model is not identifiable when multiple combinations of parameter values are solutions to the equation we are interested in solving \hyperref[csl:58]{(Ogle \& Barber, 2020}; \hyperref[csl:59]{Wieland et al., 2021)}. Identifiability issues in some complex SDMs cannot be proven mathematically but through numerical integration of its Fisher information matrix \hyperref[csl:44]{(Dorazio, 2014)}. In this study, we also check for identifiability issues in the proposed framework by describing our framework as a generalisation of the developed SDMs that account for sources of biases. Due to these identifiability issues, the inference is non-trivial and thus, the inference is left to a future paper.

This paper is organised as follows: in section 2, we introduce the multispecies biodiversity data collection process. In section 3, we explore the properties of the proposed model for biodiversity data by exploring potential issues of identifiability that this model can face and discuss possible ways to solve these issues. Finally, in section 4, we make concluding remarks and suggest future work.

\section*{Methods}

\subsection*{Motivating Example} \label{motivating}
We assume for our process model that in our study region (defined by mountains, roads, parks, and forests) with a closed system (no immigration and emigration), there are two species, each with an average number of observations (mean intensity) $\lambda_{i, true}$ for $i=1,2$ species. Some environmental covariates can explain this expected number of observations. We assume that species $1$ is more abundant (i.e. $\lambda_{i, true}$ is larger) in parks, forests and roadsides than in mountainous areas. In contrast, the opposite is true for species $2$, as presented in Figure \ref{fig:motivating}a).

Observers are likely to collect data at the most accessible places, such as roads and parks, leaving the less accessible sites, such as mountains and forests, with less data collected. We assume that the observer would collect data on both species at the same sampling location and thus the same probability of sampling, $\zeta(\textbf{s} )$, for both species, as presented in Figure \ref{fig:motivating}b).

Once observers have visited a site, there is a chance that some species are observed. This brings about imperfect detection, reducing the species distribution at the sampled locations. The probability of detecting an individual species at location $s$, $\psi(\textbf{s} )$, differs for every sampling location. In this example, we assume that the detection of species $1$ and $2$ is higher at the parks and roadside than the mountains and forests, as shown in Figure \ref{fig:motivating}c).

The detected individuals are then reported to the biodiversity database but are subject to misclassifications. Since we have a closed system, misclassification does not reduce the total number of observed individuals. However, it alters the observed amount of individuals of each species reported (i.e. the observed abundance of each species can either reduce or increase). This explanation is presented in Figure \ref{fig:motivating}d) with an equal chance of being classified as either species $1$ or $2$ when the true species identity is species $2$ and no probability of misclassification when the true identity of species is Species $1$.

The observed intensity of species $1$ and $2$, $\lambda_{obs}$, are very low and almost the same in mountainous areas. However, species $2$ has a high true intensity in that region (Figure \ref{fig:motivating}d)). Species $1$ has a high intensity in the forest but has approximately the same observed intensity as species $2$. The observed intensity is lower along the road and in the parks than the true species intensity. This explanation shows how the observed intensity of each species distribution is a thinning of the actual species distribution by sampling, detection and classification processes.\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/generatingProcess/generatingProcess}
\caption{{\label{fig:motivating} 
    (a) The true intensity $\lambda_{true}(\textbf{s})$, (b) sampling probability $\zeta(\textbf{s})$, (c) detection probability $\psi(\textbf{s})$ and (d) observed intensity $\lambda_{obs}(\textbf{s})$  in four regions (x-axis): mountains, parks, forests and roadside for two species: species 1 (solid red lines) and species 2 (dashed green lines). In (d), the true intensity for each species is coloured in grey for comparison with the observed intensity of the species. For further details of this notation, see Table \ref{tab:parameterDefintions}.
    %The observed intensity of two species after thinning the true species intensity with the probability of sampling, detection and misclassification in four regions used as x-label axis: mountains, parks, forests and roadside. The true intensity for each species $\lambda_{true}$(\textbf{s}) is plotted together with the observed intensity $\lambda_{obs}$(\textbf{s}) in d) (the true intensity is coloured in grey).%
}}
\end{center}
\end{figure}

\subsection*{Definition of terminology}
Table \ref{tab:parameterDefintions} presents the notation and definition of the parameters used for the rest of this paper.\selectlanguage{english}
\begin{table}[htbp]
%
%\tablestyle[sansbold]
%
\begin{tabular}{|b{0.2\linewidth}|b{0.7\linewidth}|}
%\theadstart
   Notation &  Definition \\

%\tbody
 %
        $\lambda_{i,obs}(\textbf{s})$ & Observed intensity (expected number of observations) of the species $i$ at location $\textbf{s} \in D \subset \mathbb{R}^2$ \\ \hline
        $\lambda_{i, true}(\textbf{s})$ &  True intensity (true expected number) of species $i$ at location $\textbf{s}$.\\ \hline

          $\zeta(\textbf{s})$ & Probability of sampling at location \textbf{s} by an observer (assumed to be the same for all the species in this study) \\ \hline

          $\psi_i(\textbf{s})$ & Probability of detection of species $i$ at location $\textbf{s}$ given the observer sampled it during the sampling process. \\ \hline

           $\Omega_{ji}$ & Probability of reporting species $i$ as species $j$ (the species is correctly classified if $i = j$). The whole matrix $\boldsymbol \Omega$ defines the classification probabilities \\ \hline

           $N$ & Number of species \\ \hline

           $\boldsymbol{\beta_{i}}$ & Parameters of the true ecological state for species $i$ \\ \hline

           $\boldsymbol{\alpha}$ & Parameters of the sampling process model \\ \hline

           $\boldsymbol{\gamma_{i}}$ & Parameters of the detection process for species $i$ \\ \hline
 %
 %\tsubheadstart
 %\tsubhead subhead &
% \tsubhead subhead \\
 %
 %content  & content \\
 %content  & content \\
 %\tend
\end{tabular}
\caption{{Definitions of the various parameters used in this paper.}}
\label{tab:parameterDefintions}

\end{table}

\subsection*{Model framework for biodiversity data}

We now develop a Bayesian framework for how the data was generated in the motivating example in section \ref{motivating}. We propose that observed biodiversity data is obtained from a thinning process on the true point pattern of species occurrences. This generating process for biodiversity data involves more than one thinning stage, with each thinning stage being a source of bias in the biodiversity data. The various stages of thinning are discussed in the following subsections and summarised in Figure \ref{fig:flowchart}.\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/flowchart/flowchart}
\caption{{\label{fig:flowchart} Flowchart of the generating process of biodiversity data described in this paper. The generating process assumes that the sampling probability, detection probability and classification probability thins the true ecological state of the species. The arrows indicate the thinning process and the dependency on the Bayesian framework used for the modelling framework. Images used in (1) - (3) were obtained from open-source sites:  \url{https://unsplash.com/s/photos/park}; \url{https://unsplash.com/s/photos/park}; \url{https://unsplash.com/s/photos/park} respectively. The image used for (4) is a screenshot of our attempted report of observation through the artsorakel app \protect\hyperref[csl:60]{(Centre, 2022)} and (5) is a screenshot of the Artsdatabanken database (\url{https://www.artsdatabanken.no/})%
}}
\end{center}
\end{figure}

As can be obtained from biodiversity databases, the biodiversity observations can be represented in the modelling framework defined by the generating process. Having a conceptual understanding of this generating process can help in modelling all potential sources of biases that can be perceived in the generating process. A graphical representation of each of the thinning stages in the flowchart is displayed in Figure \ref{csdatagpfig}.

The observed intensity can be described mathematically as:

\begin{equation}
\begin{align}
  \text{Observed intensity of species j}  &= \sum_{i=1}^{\text{N}} \text{True intensity of species i} \\
  &\times \text{P(individual i was sampled given its true intensity)} \\
  & \times \text{P(species i was detected given it was sampled)} 
  \\
  & \times \text{P(species j was reported given the true species was species i)}.
  \\
     \text {and using definitions in Table \ref{tab:parameterDefintions}:} \\ \lambda_{j,obs}(\textbf{s} ) &= \sum_{i=1}^{\text{N}} {\lambda_{i, true}(\textbf{s} )} \times {\zeta(\textbf{s} )} \times {\psi_i(\textbf{s} )} \times \Omega_{ji}.
         
\end{align}
\label{overall1} 
\end{equation}

In matrix notation, the model in equation \eqref{overall1} can be represented as:

\begin{equation}

\begin{align}

\begin{bmatrix}
          \lambda_{1,obs}(\textbf{s} )  \\
          \lambda_{2,obs}(\textbf{s} )  \\
          \vdots\\
          \lambda_{N,obs}(\textbf{s} )\\
     \end{bmatrix}
 &=

\begin{bmatrix}
          \Omega_{11}& \Omega_{21} & \ldots & \Omega_{N,1} \\
          \Omega_{12}& \Omega_{22} & \ldots & \Omega_{N,2} \\
          \vdots & \vdots & \ddots & \vdots \\
          \Omega_{1,N}& \Omega_{2,N} & \ldots & \Omega_{N,N} \\
     \end{bmatrix}
^T

\begin{bmatrix}
          \lambda_{true,1}(\textbf{s} ) \times \zeta(\textbf{s} ) \times \psi_{1}(\textbf{s} )   \\
          \lambda_{true,2}(\textbf{s} ) \times \zeta(\textbf{s} ) \times \psi_{2}(\textbf{s} ) \\
          \vdots\\
          \lambda_{true,N}(\textbf{s} ) \times \zeta(\textbf{s} ) \times \psi_{N}(\textbf{s} ),\\
     
\end{bmatrix}

\end{align}
\label{overall2} 

\end{equation}
where the parameters are defined in Table \ref{tab:parameterDefintions}.\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/csdgpfig/csdgpfig}
\caption{{\label{csdatagpfig} Illustrative example of the model framework described in section (3.3). The Figure shows the reported occurrences for each species can be downloaded from a biodiversity database. Each row represents a single species structural data collection process defined by each column.%
}}
\end{center}
\end{figure}

From equations (1) and (2), the true ecological state is affected by the sampling process, the detection probability as well as the classification probabilities and this results in species-specific point patterns with intensities $\lambda_{i,obs}, i=\{1,\ldots,N\}$. In equation \ref{overall2}, it is noteworthy that $\lambda_{i,obs}$ are weighted sums of the observed intensities before misclassification ($\lambda_{true,j}(\textbf{s} ) \times \zeta(\textbf{s} ) \times \psi_j(\textbf{s} )$) with weights given by the classification probabilities for species $i$ ($\Omega_{.i}=[\Omega_{i1},\ldots,\Omega_{i,N}]$). Hence, misclassification can result in both fewer or more reported individuals.

\subsubsection*{True ecological state process}
We assume that the true ecological state can either be modelled in a continuous space as a continuous surface that represents the expected species abundance (formally termed intensity) or in the discrete space as the expected species abundance at a site.

In this study, we model the intensity of the true ecological state of each species, $\lambda_{i, true}(\textbf{s} )$, as a Log Gaussian Cox Process \hyperref[csl:61]{(Illian et al., 2008)}:

\begin{equation}
\label{EcoLGCP} 
   \log(\lambda_{i,true}(\textbf{s} ))= X(\textbf{s} )\boldsymbol{\beta_{i}} +  \omega_{1i}(\textbf{s} ) \qquad \text{for  $i=\{1,\ldots,N\}$},

\end{equation}
where $\omega_{1i}(\textbf{s} )$ is a zero-mean Gaussian Random Field \hyperref[csl:62]{(Banerjee et al., 2015)} that differs for each species, $\boldsymbol \beta_{i}$ are the fixed effects that describe the mean of the GRF for species $i$, $\mathbf{X(\textbf{s})}$ is the set of covariates that affects the true intensity of the species and $N$ is the number of species. Given the same study region, the covariate matrix $\mathbf{X(\textbf{s})}$ is assumed, in principle, to be the same for all the species in this presentation. We assume a different GRF for each species as we assume the spatial autocorrelation varies from species to species, as presented in the left-side panel of Figure \ref{csdatagpfig}.

\subsubsection*{Uneven sampling effort}
Biodiversity data is affected by where observers prefer to go and collect observations, as discussed in the motivation example. This preference in sampling introduces a bias and is regarded as a thinning source of the ecological process. The thinning process can be modelled in the continuous space as a continuous surface that represents the probability of retaining an observation in the true ecological state or in the discrete space as the probability of collecting data at a site. The thinning source can be made to depend on the true state, for example, when the presence of the species influences where observers go to collect observation \hyperref[csl:32]{(Diggle et al., 2010)}.

In this framework, we assume that this thinning source depends on another GRF, $\omega_2(\textbf{s} )$, which we assume is the same for all the species. We do this because we assume the sampling process describes the probability that an individual is located, irrespective of what species it is. In scenarios where the target species influences the sampling process (such as when the intensity for one species influences sampling), species-specific GRFs can be defined.

The probability of retaining a point in the true ecological state (that is, an observer chooses the location occupied by an individual species) is given by:

\begin{equation}
\label{SampGLM} 
    logit(\zeta(\textbf{s} )) = U(\textbf{s} )\boldsymbol{\alpha} + \omega_2(\textbf{s} )

\end{equation}
where $U(\textbf{s} )$ are the covariates explaining the sampling process of citizen scientists, $\boldsymbol \alpha$ are the parameters associated with these covariates and $\omega_2(\textbf{s} )$ is a zero-mean GRF. An example of this thinning stage can be seen in the second column of Figure \ref{csdatagpfig}. Most of the true occurrences in the blue area are not reported for both species. Hence, the number of occurrences reported after the first stage of thinning is about 61\% of the true abundance for species 1 and 73\% for species 2.

\subsubsection*{Detection process}
Given the true ecological state and the sampling process, an individual species may be present at a given sampled location. Still, they might not be observed, reducing the observed species abundance. Hence, this model framework incorporates the detection probability by assuming a generalised linear model for the detection. The detection probability model is given as:

\begin{equation}
\label{DetectGLM} 
    logit(\psi_{i}(\textbf{s} )) = Z(\textbf{s} )\boldsymbol{\gamma_i},

\end{equation}
where $Z(\textbf{s} )$ are the covariates that affect the detection of the species, such as the weather, size of species, and distance from the observer to observation, and $\gamma$ are the fixed effects that describe the covariate effects for each species. In the third column of Figure \ref{csdatagpfig}, we can see how the detectability of each species affects the resulting point pattern. For example, species $2$ has high detectability almost everywhere in our study region. In contrast, species $1$ has low detectability in the upper right corner, which means the loss of many occurrences because they could not be detected. Consequently, the number of occurrences after the second thinning stage is about 31\% of the sampled species $1$ occurrences and 23\% of the sampled species $2$.

\subsubsection*{Classification process}

Observers report the observations they have detected to a database they have signed up to. These detected observations are likely to be misidentified as another species or misreported to another species, as have been noted and modelled in some SDMs \hyperref[csl:47]{(Adjei et al., 2023}; \hyperref[csl:50]{Wright et al., 2020}; \hyperref[csl:49]{Chambert et al., 2015)}. By this, we assume that observers are sure to report all their detections. These misidentifications arise due to factors such as the observer's experience and easy confusion about similar-looking species while the misreporting could arise from the observer's difficulty in reporting species to the database. We collectively refer to both misidentification and misreporting as misclassifications in this framework. We account for these misclassifications in our model framework as:

\begin{equation}
Observation_i \sim Multinomial (\Omega_{i.} )

\end{equation}
where $\Omega_{ij} =  P(\text{classifying species j}| \text{true species i})$ is the classification probabilities (species is correctly classified when $i = j$). Given the definition of $\Omega$, the rows of this matrix should sum to 1. The classification probabilities can be modelled as constant in space and drawn from a Dirichlet distribution \hyperref[csl:50]{(Wright et al., 2020}; \hyperref[csl:63]{Chambert et al., 2018}; \hyperref[csl:64]{Chambert et al., 2018)} or heterogeneous across the space by either fitting a multinomial generalised linear model to each of the rows of $\Omega$ or using ML prediction scores as a proxy to the classification probabilities \hyperref[csl:47]{(Adjei et al., 2023)}. The last column of Figure \ref{csdatagpfig} shows the reported species with the potential misclassifications that observers can make. The misclassification does not reduce the total observed species; however, the number of observed occurrences for species $1$ reduced by $5$ and that of species $2$ increased by $5$.

\section*{Identifiability issues}

The model framework we have defined in section $3.3$ may be overparameterized (i.e. have more model parameters than can be estimated from the data). Overparameterized models are most likely to have identifiability issues \hyperref[csl:59]{(Wieland et al., 2021}; \hyperref[csl:65]{Little et al., 2010)} because a change in one of the parameters compensates exactly for the changes in the posterior distribution caused by a change in another parameter \hyperref[csl:66]{(Rannala, 2002}; \hyperref[csl:67]{Swartz et al., 2004}; \hyperref[csl:68]{Raue et al., 2013}; \hyperref[csl:58]{Ogle \& Barber, 2020)}. In the proposed framework, non-identifiability refers to different combinations of parameters producing the same posterior probability or likelihood \hyperref[csl:66]{(Rannala, 2002}; \hyperref[csl:68]{Raue et al., 2013}; \hyperref[csl:59]{Wieland et al., 2021)} or the inability of models, data or implementation to inform about effects of interest \hyperref[csl:68]{(Raue et al., 2013}; \hyperref[csl:58]{Ogle \& Barber, 2020)}. These identifiability issues can lead to poor mixing and convergence issues when Markov chain Monte Carlo (MCMC) methods are used \hyperref[csl:58]{(Ogle \& Barber, 2020)} since the MCMC algorithms would get stuck trying to explore a ridged likelihood \hyperref[csl:67]{(Swartz et al., 2004)}. This results in correlated and/or biased estimated parameters, and thus biased inferences of species distributions \hyperref[csl:67]{(Swartz et al., 2004}; \hyperref[csl:58]{Ogle \& Barber, 2020)}. Because of the challenges posed by these identifiability issues, we assessed the potential sources of non-identifiability in our proposed modelling framework.

In this section, we provide intuitive scenarios when identifiability issues can occur in this model, using our motivation example (see Section 3.1) as the baseline scenario (summaries of the scenarios are presented in Table \ref{tab:identifiability}). The structural identifiability happens in our framework model when, for example, there are no misclassifications and the true intensity, $\lambda_{true}(\textbf{s} )$ is multiplied by a constant $K \in \mathbb{R}$, and either the sampling retaining probabilities $\zeta(\textbf{s} )$ or the detection probabilities $\psi_i(\textbf{s} )$ are multiplied by a constant $1/K$. That is, either

\begin{equation}

\begin{align}
     \lambda_{i,obs} &= K \cdot {\lambda_{i, true}(\textbf{s} )} \times \bigg(\frac{1}{K}\cdot {\zeta(\textbf{s} )} \bigg) \times {\psi_i(\textbf{s} )}  \quad \text{or}\\
     \lambda_{i,obs} &= K \cdot {\lambda_{i, true}(\textbf{s} )} \times  {\zeta(\textbf{s} )} \times \bigg(\frac{1}{K}\cdot {\psi_i(\textbf{s} )}\bigg) 
     
\end{align}
     \label{identifiability} 

\end{equation}
would produce the same observed intensity $\lambda_{i,obs}$.

Assuming perfect classifications (ie. $\Omega_{12} =\Omega_{21}=0$), it can be observed that scaling down the true intensity of species $1$ ($\lambda_{1, true}$) by a factor of $0.5$ and scaling up the probability of sampling along the four locations ($\zeta(\textbf{s} )$) by a factor of $0.5$ as well as multiplying the true intensity of species $1$ ($\lambda_{1, true}$) by a factor of $0.5$ and dividing its detection probability ($\psi_{1}(\textbf{s} )$) at all four locations by $0.5$ gives the same observed intensity ($\lambda_{1, obs}$) as the baseline scenario. When there are misclassifications ($\Omega_{12}$ and $\Omega_{21}$ are different from zero), then scaling down the true intensity of species $1$ ($\lambda_{1, true}(\textbf{s} )$) by a factor of $0.5$ and dividing its detection probability ($\psi_{1}(\textbf{s} )$) by $0.5$ as well as multiplying the true intensity of both species $1$ ($\lambda_{1, true}(\textbf{s} )$) and $2$ ($\lambda_{2, true}(\textbf{s} )$) by a factor of $0.5$ and dividing the sampling probability ($\zeta (\textbf{s} )$) by $0.5$ gives the same observed intensity as the baseline scenario for species $1$ ($\lambda_{1, obs}(\textbf{s} )$) in both instances and ($\lambda_{2, obs}(\textbf{s} )$) for the latter.\selectlanguage{english}
\begin{table}
  \centering

\begin{tabular}
{cc|ccccccccc|cc}
\cmidrule{3-13}          &       & $\lambda_{1,true}(s)$ & $\lambda_{2,true}(s)$ & $\zeta(s)$ & $\psi_1(s)$ & $\psi_2(s)$ & $\Omega_{11}$ & $\Omega_{12}$ & $\Omega_{21}$ & $\Omega_{22}$ & $\lambda_{1,obs}(s)$ & $\lambda_{2,obs}(s)$ \\
    \midrule
    \multirow{3}[2]{*}{No misclassification} & Baseline scenario & 10    & 7     & 0.5   & 0.2   & 0.5   & 1     & 0     & 0     & 1     & \textcolor[rgb]{ 1,  0,  0}{1} & 1.75 \\
          & $K \times \lambda_{1,true}$ and $\zeta/K$ & 5  & 7     & 1 & 0.2   & 0.5   & 1     & 0     & 0     & 1     & \textcolor[rgb]{ 1,  0,  0}{1} & 0.0175 \\
          & $K \times \lambda_{1,true}$ and $\psi_1/K$ & 5  & 7     & 0.5   & 0.4 & 0.5   & 1     & 0     & 0     & 1     & \textcolor[rgb]{ 1,  0,  0}{1} & 1.75 \\
    \midrule
    \multirow{3}[2]{*}{Misclassification} & Baseline scenario & 10    & 7     & 0.5   & 0.2   & 0.5   & 0.5   & 0.5   & 0.2   & 0.8   & \textcolor[rgb]{ 1,  0,  0}{0.85} & \textcolor[rgb]{ 1,  0,  0}{1.9} \\
          & $K \times \lambda_{1,true}$ and $\psi_1/K$ & 5  & 7     & 0.5   & 0.4& 0.5   & 0.5   & 0.5   & 0.2   & 0.8   & \textcolor[rgb]{ 1,  0,  0}{0.85} & 1.9 \\
          & $K \times \lambda_{1,true}$, $K \times \lambda_{2,true}$ and $\zeta/K$ & 5  & 3.5  & 1 & 0.2   & 0.5   & 0.5   & 0.5   & 0.2   & 0.8   & \textcolor[rgb]{ 1,  0,  0}{0.85} & \textcolor[rgb]{ 1,  0,  0}{1.9} \\
    \bottomrule
    \end{tabular}
      \caption{{Structural identifiability scenarios with and without misclassification for two species. In red: Equal observed intensities $\lambda_{obs}$ with different parameter values with $K=0.5$.}}
  \label{tab:identifiability} 
\end{table}

Although our intuitive examples have shown identifiability issues with changes in true ecological state, sampling and detection process parameters, species distribution models that estimate classification probabilities together with ecological parameters have noted non-identifiability of parameters if there is no additional information available to uniquely estimate the classification probabilities \hyperref[csl:69]{(Ferguson et al., 2015}; \hyperref[csl:70]{Guillera-Arroita et al., 2017}; \hyperref[csl:71]{Clare et al., 2021)}.

Notably, these challenges associated with identifiability can be revealed through the MCMC chains after fitting and estimating the parameters in the model \hyperref[csl:72]{(Gelman \& Hill, 2006}; \hyperref[csl:58]{Ogle \& Barber, 2020)}. The fitting and estimation of the parameters are out of the scope of this paper. It is possible that all the model parameters may not be identifiable, which would not be an issue if inferences on those parameters are not interesting. However, we need to be aware of these problems and why it is important, how to diagnose them and the possible ways to handle them. We propose these potential identifiability issues can be mitigated by integrating different data sources for the various parts of the model framework or using informed priors if multiple datasets from different sampling protocols cannot be obtained. These mitigating methods are discussed in the subsections that follow.

\subsection*{Informed priors}

More generally, the introduction of informed priors on model parameters can help overcome identifiability issues \hyperref[csl:73]{(Gelman et al., 1996}; \hyperref[csl:74]{Zyphur \& Oswald, 2015}; \hyperref[csl:75]{Lemoine, 2019}; \hyperref[csl:67]{Swartz et al., 2004}; \hyperref[csl:58]{Ogle \& Barber, 2020)}. These informative priors represent our prior knowledge about the effects of biases on the true ecological state. These informed priors improve MCMC performance and estimation of the model parameters \hyperref[csl:58]{(Ogle \& Barber, 2020)}.

The informed priors can be obtained from expert opinions and elicitation techniques \hyperref[csl:76]{(Murray et al., 2009}; \hyperref[csl:77]{Martin et al., 2005}; \hyperref[csl:78]{Kuhnert, 2011)}, previously fitted models or algorithms \hyperref[csl:47]{(Adjei et al., 2023}; \hyperref[csl:79]{McDonald \& Hodgson, 2021)}, or results from the literature and historical records \hyperref[csl:80]{(Stirling et al., 2016)}.

Although proposing informed priors on model parameters can result in technically identifiable models, it can be hard to develop such prior distributions \hyperref[csl:75]{(Lemoine, 2019)}. Moreover, the results will most likely depend on the choice of the prior distribution, even with larger sample sizes \hyperref[csl:58]{(Ogle \& Barber, 2020)}, consequently leading to poor inferences. Therefore, these priors have to be chosen with care.

\subsection*{Data Integration}
 Another approach to inform the model parameters and ensure identifiability is to provide additional data sources for some or all the processes described in the framework. This approach relies on data integration, which involves integrating data from different sources and quality to estimate shared parameters \hyperref[csl:81]{(Miller et al., 2019}; \hyperref[csl:46]{Isaac et al., 2020)}.

The joint likelihood of our Bayesian framework can be expressed as a combination of likelihoods from different independent data sources. That is, we can express our likelihood as:

\begin{equation}
\label{IDMLikelihood}

\begin{align}
    L(\boldsymbol{\beta}, \boldsymbol{\alpha}, \boldsymbol{\gamma}, \boldsymbol{\Omega}) =   L(\boldsymbol{\beta}) \times  L(\boldsymbol{\alpha}|\boldsymbol{\beta}) \times    L( \boldsymbol{\gamma}|\boldsymbol{\beta})
    \times    L(\boldsymbol{\Omega}|\boldsymbol{\beta});
    
\end{align}

\end{equation}
where $L(\boldsymbol{\beta}) \times  L(\boldsymbol{\alpha}|\boldsymbol{\beta})$ is the expression of the joint likelihood of integrated distribution models (IDMs, hereafter) that account for sampling effort whilst estimating true ecological state parameters, $L(\boldsymbol{\beta}) \times  L(\boldsymbol{\gamma}|\boldsymbol{\beta})$ is the expression of the joint likelihood of IDMs that account for imperfect detection whilst estimating true ecological state parameters and $L(\boldsymbol{\beta}) \times  L(\boldsymbol{\Omega}|\boldsymbol{\beta})$ is the expression of the joint likelihood of IDMs that account for misclassification whilst estimating true ecological state parameters. There is also the possibility that a given data source can inform the estimation of two sources of biases whilst estimating the true ecological state parameters. For example, repeated-visit occupancy data collected using preferential sampling can inform the estimation of the true occupancy probability, detection probability and sampling bias \hyperref[csl:82]{(Case \& Lawler, 2017}; \hyperref[csl:44]{Dorazio, 2014)}; repeated-visit occupancy data with either site or observation confirmation can inform the estimation of the occupancy probability, detection probability and classification probabilities \hyperref[csl:69]{(Ferguson et al., 2015}; \hyperref[csl:71]{Clare et al., 2021)}. When one dataset informs the estimation of two sources of biases, the joint likelihood will also be a combination of the joint likelihood framework expressed in equation \eqref{IDMLikelihood}, whether we assume the two biases are independent of each other or not (shown by the intersection of two sources of biases in Figure \ref{fig:dataIntegration}). We assume that no single dataset can inform the estimation of all the model parameters. The possible combinations of informing the estimation of model parameters in our proposed framework are presented in Figure \ref{fig:dataIntegration}. \\\selectlanguage{english}
\begin{figure}[H]
\begin{center}
\includegraphics[width=0.70\columnwidth]{figures/dataIntegration-chat/dataIntegration-chat}
\caption{{\label{fig:dataIntegration} Venn diagram showing the possible ways of integrating data to estimate the model parameters in our proposed framework. For each dataset obtained, we assume it should be geared towards estimating the true ecological process parameters ($\boldsymbol{\beta}$), but provides information to estimate at least one of the other process parameters: sampling process ($\boldsymbol{\alpha}$), detection process ($\boldsymbol{\gamma}$)  and classification process ($\boldsymbol{\Omega}$). The intersection of two processes denotes data that can inform the estimation of the true ecological state and the two parameters in question.%
}}
\end{center}
\end{figure}

We present potential data sources that can be integrated into this biodiversity data collection framework by gaining insight from some IDMs that have accounted for the sources of bias discussed in this paper.  IDMs that account for spatial biases or those that use data obtained from preferential sampling protocols \hyperref[csl:83]{(Gelfand \& Shirota, 2019)} have integrated presence-only data from opportunistic surveys with either presence-absence data \hyperref[csl:30]{(Simmonds et al., 2020}; \hyperref[csl:84]{Fletcher et al., 2016}; \hyperref[csl:43]{Koshkina et al., 2017)}, detection/non-detection data \hyperref[csl:85]{(Adde et al., 2021)} or camera trap and telemetry data \hyperref[csl:86]{(Cretois et al., 2021}; \hyperref[csl:27]{Sicacha-Parada et al., 2021}; \hyperref[csl:87]{Santangeli et al., 2020)} to gain insight about the sampling process whilst making inferences about the ecological state of interest. IDMs that account for imperfect detection have used data from repeated counts \hyperref[csl:44]{(Dorazio, 2014}; \hyperref[csl:10]{K{\'e}ry \& Royle, 2015}; \hyperref[csl:11]{K{\'e}ry \& Royle, 2020)}, double observers or site occupancy protocols \hyperref[csl:88]{(Canessa et al., 2012}; \hyperref[csl:43]{Koshkina et al., 2017)} to inform the estimation of the model parameters. The repeated surveys in the data provide information about the abundance and detection process, making the parameters identifiable \hyperref[csl:44]{(Dorazio, 2014)}. IDMs that account for misclassification need extra information from the verification process to make the parameters identifiable. These verifications can be done on the site-level \hyperref[csl:63]{(Chambert et al., 2018)}, aggregated individual sample level  \hyperref[csl:50]{(Wright et al., 2020)}, or individual sample-level \hyperref[csl:89]{(Spiers et al., 2022)}. It must be stated that some studies have also accounted for multiple biases such as misclassification and detection probability \hyperref[csl:69]{(Ferguson et al., 2015}; \hyperref[csl:90]{Louvrier et al., 2018}; \hyperref[csl:91]{Clement et al., 2022}; \hyperref[csl:42]{Royle \& Link, 2006}; \hyperref[csl:48]{Miller et al., 2011)}, sampling bias and detection probability \hyperref[csl:44]{(Dorazio, 2014)}.

With these potential identifiability issues in our model, additional data is needed to constrain the parameters of each thinning stage in our model. Integrating the various sources of data described in the previous paragraph takes advantage of the insights that professional repeated surveys provide about the detection process, a distinct covariate such as distance to the road about the sampling process, presence-absence data regarding sampling bias and data on the verification process for the misclassification probabilities.

\section*{Conclusion}

Biodiversity data have many data quality issues \hyperref[csl:92]{(Dickinson et al., 2012}; \hyperref[csl:7]{Bal{\'a}zs et al., 2021)}. These issues include but are not limited to sampling bias, imperfect detection and misclassification. When these biases are left unaccounted for, they can affect the inferences made about the species distribution. In this paper, we have proposed a multi-species approach that makes it possible to account for these common sources of biases in biodiversity data.  This general model framework assumes these sources of biases act as a thinning process on the true ecological state.

By modelling the sources of biases in biodiversity data, we also described a structural model for the process of collecting biodiversity data. Existing frameworks have accounted for some sources of biases in their SDMs (as we have discussed in section (4) of this paper). We add to this knowledge by providing a general hierarchical structure that describes how sources of bias in biodiversity data can be modelled probabilistically. For example, other sources of biases, such as taxonomic bias, and geographic biases, can be modelled with this framework. Since failure to properly account for biases in species distributions can lead to poor inferences of species distributions \hyperref[csl:44]{(Dorazio, 2014}; \hyperref[csl:93]{Royle \& Nichols, 2003}; \hyperref[csl:11]{K{\'e}ry \& Royle, 2020)}, this novel approach will allow scientists to better understand challenges in biodiversity data and how to handle them.

In this study, we only accounted for uneven sampling effort, imperfect detections and misclassifications as biases in biodiversity data. However, we acknowledge that other relevant sources of biases may be explicitly ignored in this framework. For example, we only described spatial biases in the framework, but temporal biases  \hyperref[csl:94]{(Courter et al., 2013)}, taxonomic biases \hyperref[csl:95]{(Di Cecco et al., 2021)}, reporting biases \hyperref[csl:57]{(Arazy \& Malkinson, 2021)} etc may be of interest. These sources of biases can be included in the proposed framework by carefully considering how they thin the true ecological process and the information needed to estimate the parameters of that model. For example, metadata from experience-based protocols of collecting biodiversity data (where the expertise of the observer can affect species identification) can help in estimating reporting probabilities.

Structural identifiability issues arise naturally in our modelling framework as the observed intensity depends on many parameters and multiple parameter combinations give the same observed intensity \hyperref[csl:59]{(Wieland et al., 2021)}. There are several instances where these identifiability issues can arise such as having correlated covariates and intercept-only models for the detection and sampling process. We have proposed in this study that these identifiability issues can be succumbed by integrating data from different sources or using informed priors. Biodiversity data collected from different sampling protocols such as opportunistic, distance sampling, and repeated count surveys can be modelled together, and each dataset provides information in this IDM to reduce the bias in parameter estimates \hyperref[csl:43]{(Koshkina et al., 2017)} or weaknesses of each dataset \hyperref[csl:54]{(Isaac et al., 2014}; \hyperref[csl:84]{Fletcher et al., 2016)}.

Our framework can be used to develop models that can be fitted to data; this will be explored in a later paper. The flexibility of the framework means that it can be used across datasets, and will help isolate the different biases that affect the observation process so they can be efficiently analysed and understood.

\section*{Acknowledgements}
This work is part of the Transforming Citizen Science for Biodiversity project, funded by the
NTNU digital transformation initiative .

\section*{Data accessibility}

The Rcode used in generating the plots in the paper can be accessed openly on \hyperref[csl:96]{(Adjei et al., 2022)} and GitHub repository \url{https://github.com/Peprah94/Data-generating-process}.

\section*{Author contributions}

All authors were involved in the framework development. KPA and JSP led the writing of the manuscript, and all authors contributed to writing the manuscript.

\selectlanguage{english}
\FloatBarrier
\section*{References}\sloppy
\phantomsection
\label{csl:1}{Statistical solutions for error and bias in global citizen science datasets}. (2014). \textit{Biological Conservation}, \textit{173}, 144–154.

\phantomsection
\label{csl:2}{iNaturalist}. (2014). \textit{Reference Reviews}.

\phantomsection
\label{csl:3}{Biodiversity information goes public: GBIF at your service}. (2011). \textit{Nordic Journal of Botany}, \textit{29}(3), 378–381.

\phantomsection
\label{csl:4}{The eBird enterprise: An integrated approach to development and application of citizen science}. (2014). \textit{Biological Conservation}, \textit{169}, 31–40.

\phantomsection
\label{csl:5}{Influence of volunteer and project characteristics on data quality of biological surveys}. (2015). \textit{Conservation Biology}, \textit{29}(3), 713–723.

\phantomsection
\label{csl:6}{Assessing citizen science data quality: an invasive species case study}. (2011). \textit{Conservation Letters}, \textit{4}(6), 433–442.

\phantomsection
\label{csl:7}{Data quality in citizen science}. (2021). \textit{The Science of Citizen Science}, 139.

\phantomsection
\label{csl:8}{Marine invasive species: validation of citizen science and implications for national monitoring networks}. (2008). \textit{Biological Invasions}, \textit{10}(1), 117–128.

\phantomsection
\label{csl:9}{Assessing the changing flowering date of the common lilac in North America: a random coefficient model approach}. (2012). \textit{Geoinformatica}, \textit{16}(4), 675–690.

\phantomsection
\label{csl:10}\textit{{Applied Hierarchical Modeling in Ecology: Analysis of Distribution, Abundance and Species Richness in R and BUGS}}. (2015).

\phantomsection
\label{csl:11}\textit{{Applied Hierarchical Modeling in Ecology: Analysis of distribution, abundance and species richness in R and BUGS: Volume 2: Dynamic and Advanced Models}}. (2020). Academic Press.

\phantomsection
\label{csl:12}{Evaluating Model Fit and Selecting among Multiple Models}. (2020). In \textit{Joint Species Distribution Modelling: With Applications in R} (pp. 217–252). Cambridge University Press. \url{https://doi.org/10.1017/9781108591720.012}

\phantomsection
\label{csl:13}{Maximum entropy modeling of species geographic distributions}. (2006). \textit{Ecological Modelling}, \textit{190}(3), 231–259. \url{https://doi.org/https://doi.org/10.1016/j.ecolmodel.2005.03.026}

\phantomsection
\label{csl:14}{Sample selection bias and presence-only distribution models: implications for background and pseudo-absence data}. (2009). \textit{Ecological Applications}, \textit{19}(1), 181–197. \url{https://doi.org/10.1890/07-2153.1}

\phantomsection
\label{csl:15}{Point pattern modelling for degraded presence-only data over large regions}. (2011). \textit{Journal of the Royal Statistical Society: Series C (Applied Statistics)}, \textit{60}(5), 757–776.

\phantomsection
\label{csl:16}{Preferential sampling for presence/absence data and for fusion of presence/absence data with presence-only data}. (2019). \textit{Ecological Monographs}, \textit{89}(3), e01372. \url{https://doi.org/10.1002/ecm.1372}

\phantomsection
\label{csl:17}{Point process models for presence-only analysis}. (2015). \textit{Methods in Ecology and Evolution}, \textit{6}(4), 366–379.

\phantomsection
\label{csl:18}{Extended statistical approaches to modelling spatial pattern in biodiversity in northeast New South Wales. II. Community-level modelling}. (2002). \textit{Biodiversity {\&} Conservation}, \textit{11}(12), 2309–2338. \url{https://doi.org/10.1023/A:1021374009951}

\phantomsection
\label{csl:19}{Selecting pseudo-absences for species distribution models: how, where and how many?}. (2012). \textit{Methods in Ecology and Evolution}, \textit{3}(2), 327–338. \url{https://doi.org/10.1111/j.2041-210X.2011.00172.x}

\phantomsection
\label{csl:20}{Sample selection bias and presence-only distribution models: implications for background and pseudo-absence data}. (2009). \textit{Ecological Applications}, \textit{19}(1), 181–197.

\phantomsection
\label{csl:21}{The effects of sampling bias and model complexity on the predictive performance of MaxEnt species distribution models}. (2013). \textit{PloS One}, \textit{8}(2), e55158.

\phantomsection
\label{csl:22}{Predicting species spatial distributions using presence-only data: a case study of native New Zealand ferns}. (2002). \textit{Ecological Modelling}, \textit{157}(2-3), 261–280.

\phantomsection
\label{csl:23}{The ghost of past species occurrence: improving species distribution models for presence-only data}. (2006). \textit{Journal of Applied Ecology}, \textit{43}(4), 802–815.

\phantomsection
\label{csl:24}{An improved approach for predicting the distribution of rare and endangered species from occurrence and pseudo-absence data}. (2004). \textit{Journal of Applied Ecology}, \textit{41}(2), 263–274.

\phantomsection
\label{csl:25}{Assessing the effects of pseudo-absences on predictive distribution model performance}. (2008). \textit{Ecological Modelling}, \textit{210}(4), 478–486.

\phantomsection
\label{csl:26}{Bias correction in species distribution models: pooling survey and collection data for multiple species}. (2015). \textit{Methods in Ecology and Evolution}, \textit{6}(4), 424–438.

\phantomsection
\label{csl:27}{Accounting for spatial varying sampling effort due to accessibility in Citizen Science data: A case study of moose in Norway}. (2021). \textit{Spatial Statistics}, \textit{42}, 100446.

\phantomsection
\label{csl:28}{Point process models for spatio-temporal distance sampling data from a large-scale survey of blue whales}. (2017). \textit{The Annals of Applied Statistics}, \textit{11}(4), 2270–2297.

\phantomsection
\label{csl:29}{What determines spatial bias in citizen science? Exploring four recording schemes with different proficiency requirements}. (2016). \textit{Diversity and Distributions}, \textit{22}(11), 1139–1149.

\phantomsection
\label{csl:30}{Is more data always better? A simulation study of benefits and limitations of integrated distribution models}. (2020). \textit{Ecography}, \textit{43}(10), 1413–1422.

\phantomsection
\label{csl:31}{Accessibility maps as a tool to predict sampling bias in historical biodiversity occurrence records}. (2019). \textit{Ecography}, \textit{42}(1), 125–136. \url{https://doi.org/10.1111/ecog.03944}

\phantomsection
\label{csl:32}{Geostatistical inference under preferential sampling}. (2010). \textit{Journal of the Royal Statistical Society: Series C (Applied Statistics)}, \textit{59}(2), 191–232. \url{https://doi.org/10.1111/j.1467-9876.2009.00701.x}

\phantomsection
\label{csl:33}{Comparing sample bias correction methods for species distribution modeling using virtual species}. (2021). \textit{Ecosphere}, \textit{12}(3), e03422.

\phantomsection
\label{csl:34}{Imperfect observations in ecological studies}. (2016). \textit{Environmental and Ecological Statistics}, \textit{23}(3), 337–358.

\phantomsection
\label{csl:35}{Integrated species distribution models: A comparison of approaches under different data quality scenarios}. (2021). \textit{Diversity and Distributions}, \textit{27}(6), 1066–1075.

\phantomsection
\label{csl:36}{Inhomogeneous Poisson point process for species distribution modelling: relative performance of methods accounting for sampling bias and imperfect detection}. (2022). \textit{Modeling Earth Systems and Environment}, 1–14.

\phantomsection
\label{csl:37}{Monitoring programs need to take into account imperfect species detectability}. (2004). \textit{Basic and Applied Ecology}, \textit{5}(1), 65–73.

\phantomsection
\label{csl:38}{Accounting for imperfect detection in observational studies: modeling wolf sightability in Yellowstone National Park}. (2020). \textit{Ecosphere}, \textit{11}(6), e03152.

\phantomsection
\label{csl:39}{Estimating site occupancy rates when detection probabilities are less than one}. (2002). \textit{Ecology}, \textit{83}(8), 2248–2255.

\phantomsection
\label{csl:40}{A Framework of Observer-Based Biases in Citizen Science Biodiversity Monitoring: Semi-Structuring Unstructured Biodiversity Monitoring Protocols}. (2022). \textit{Citizen Science for Future Generations}.

\phantomsection
\label{csl:41}{Fitting and interpreting occupancy models}. (2013). \textit{PloS One}, \textit{8}(1), e52015.

\phantomsection
\label{csl:42}{Generalized site occupancy models allowing for false positive and false negative errors}. (2006). \textit{Ecology}, \textit{87}(4), 835–841.

\phantomsection
\label{csl:43}{Integrated species distribution models: combining presence-background data and site-occupancy data with imperfect detection}. (2017). \textit{Methods in Ecology and Evolution}, \textit{8}(4), 420–430.

\phantomsection
\label{csl:44}{Accounting for imperfect detection and survey bias in statistical analysis of presence-only data}. (2014). \textit{Global Ecology and Biogeography}, \textit{23}(12), 1472–1484.

\phantomsection
\label{csl:45}{Accounting for imperfect detection in data from museums and herbaria when modeling species distributions: combining and contrasting data-level versus model-level bias correction}. (2021). \textit{Ecography}, \textit{44}(9), 1341–1352.

\phantomsection
\label{csl:46}{Data integration for large-scale models of species distributions}. (2020). \textit{Trends in Ecology \& Evolution}, \textit{35}(1), 56–67.

\phantomsection
\label{csl:47}{Modelling heterogeneity in the classification process in multi-species distribution models can improve predictive performance}. (2023). \textit{ArXiv Preprint ArXiv:2305.01989}.

\phantomsection
\label{csl:48}{Improving occupancy estimation when two types of observational error occur: Non-detection and species misidentification}. (2011). \textit{Ecology}, \textit{92}(7), 1422–1428.

\phantomsection
\label{csl:49}{Modeling false positive detections in species occurrence data under different study designs}. (2015). \textit{Ecology}, \textit{96}(2), 332–339.

\phantomsection
\label{csl:50}{Modelling misclassification in multi-species acoustic data when estimating occupancy and relative activity}. (2020). \textit{Methods in Ecology and Evolution}, \textit{11}(1), 71–81.

\phantomsection
\label{csl:51}{Dependent double-observer method reduces false-positive errors in auditory avian survey data}. (2020). \textit{Ecological Applications}, \textit{30}(2), e02026.

\phantomsection
\label{csl:52}{Hierarchical spatial models of abundance and occurrence from imperfect survey data}. (2007). \textit{Ecological Monographs}, \textit{77}(3), 465–481.

\phantomsection
\label{csl:53}{Species identification by experts and non-experts: comparing images from field guides}. (2016). \textit{Scientific Reports}, \textit{6}(1), 1–7.

\phantomsection
\label{csl:54}{Statistics for citizen science: extracting signals of change from noisy ecological data}. (2014). \textit{Methods in Ecology and Evolution}, \textit{5}(10), 1052–1060.

\phantomsection
\label{csl:55}\textit{{Biological diversity: frontiers in measurement and assessment}}. (2010). OUP Oxford.

\phantomsection
\label{csl:56}\textit{{Do we know enough? The impacts of data biases, insufficient sampling and degradation on biodiversity estimates in Tanzanian forests--implications for conservation planning}}. (2010). [PhD thesis]. University of York.

\phantomsection
\label{csl:57}{A Framework of Observer-Based Biases in Citizen Science Biodiversity Monitoring: Semi-Structuring Unstructured Biodiversity Monitoring Protocols}. (2021). \textit{Frontiers in Ecology and Evolution}, 435.

\phantomsection
\label{csl:58}{Ensuring identifiability in hierarchical mixed effects Bayesian models}. (2020). \textit{Ecological Applications}, \textit{30}(7), e02159.

\phantomsection
\label{csl:59}{On structural and practical identifiability}. (2021). \textit{Current Opinion in Systems Biology}.

\phantomsection
\label{csl:60}\textit{{Artsorakel (Version 2.01.04) [Mobile app] Google Play.}}. (2022). Version. \url{https://play.google.com/store/apps/details?id=no.artsdatabanken.orakel}

\phantomsection
\label{csl:61}\textit{{Statistical analysis and modelling of spatial point patterns}} (Vol. 70). (2008). (Vol.). John Wiley \& Sons.

\phantomsection
\label{csl:62}\textit{{Hierarchical Modeling and Analysis for Spatial Data, Second Edition}} (2ed. ed.). (2015). (ed.). CRC Press. \url{http://gen.lib.rus.ec/book/index.php?md5=ee60a1077624ddac0ebb2042fdf2098e}

\phantomsection
\label{csl:63}{A new framework for analysing automated acoustic species detection data: Occupancy estimation and optimization of recordings post-processing}. (2018). \textit{Methods in Ecology and Evolution}, \textit{9}(3), 560–570.

\phantomsection
\label{csl:64}{Two-species occupancy modelling accounting for species misidentification and non-detection}. (2018). \textit{Methods in Ecology and Evolution}, \textit{9}(6), 1468–1477.

\phantomsection
\label{csl:65}{Parameter identifiability and redundancy: theoretical considerations}. (2010). \textit{PloS One}, \textit{5}(1), e8915.

\phantomsection
\label{csl:66}{Identifiability of parameters in MCMC Bayesian inference of phylogeny}. (2002). \textit{Systematic Biology}, \textit{51}(5), 754–760.

\phantomsection
\label{csl:67}{Bayesian identifiability and misclassification in multinomial data}. (2004). \textit{Canadian Journal of Statistics}, \textit{32}(3), 285–302.

\phantomsection
\label{csl:68}{Joining forces of Bayesian and frequentist methodology: a study for inference in the presence of non-identifiability}. (2013). \textit{Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences}, \textit{371}(1984), 20110544.

\phantomsection
\label{csl:69}{Occupancy models for data with false positive and false negative errors and heterogeneity across sites and surveys}. (2015). \textit{Methods in Ecology and Evolution}, \textit{6}(12), 1395–1406.

\phantomsection
\label{csl:70}{Dealing with false-positive and false-negative errors about species occurrence at multiple levels}. (2017). \textit{Methods in Ecology and Evolution}, \textit{8}(9), 1081–1091.

\phantomsection
\label{csl:71}{Generalized model-based solutions to false-positive error in species detection/nondetection data}. (2021). \textit{Ecology}, \textit{102}(2), e03241.

\phantomsection
\label{csl:72}\textit{{Data analysis using regression and multilevel/hierarchical models}}. (2006). Cambridge university press.

\phantomsection
\label{csl:73}{Physiological pharmacokinetic analysis using population modeling and informative prior distributions}. (1996). \textit{Journal of the American Statistical Association}, \textit{91}(436), 1400–1412.

\phantomsection
\label{csl:74}{Bayesian estimation and inference: A user’s guide}. (2015). \textit{Journal of Management}, \textit{41}(2), 390–420.

\phantomsection
\label{csl:75}{Moving beyond noninformative priors: why and how to choose weakly informative priors in Bayesian analyses}. (2019). \textit{Oikos}, \textit{128}(7), 912–928.

\phantomsection
\label{csl:76}{How useful is expert opinion for predicting the distribution of a species within and beyond the region of expertise? A case study using brush-tailed rock-wallabies Petrogale penicillata}. (2009). \textit{Journal of Applied Ecology}, \textit{46}(4), 842–851.

\phantomsection
\label{csl:77}{The power of expert opinion in ecological models using Bayesian methods: impact of grazing on birds}. (2005). \textit{Ecological Applications}, \textit{15}(1), 266–280.

\phantomsection
\label{csl:78}{Four case studies in using expert opinion to inform priors}. (2011). \textit{Environmetrics}, \textit{22}(5), 662–674.

\phantomsection
\label{csl:79}{Counting Cats: The integration of expert and citizen science data for unbiased inference of population abundance}. (2021). \textit{Ecology and Evolution}, \textit{11}(9), 4325–4338.

\phantomsection
\label{csl:80}{Using verified species distribution models to inform the conservation of a rare marine species}. (2016). \textit{Diversity and Distributions}, \textit{22}(7), 808–822.

\phantomsection
\label{csl:81}{The recent past and promising future for data integration methods to estimate species’ distributions}. (2019). \textit{Methods in Ecology and Evolution}, \textit{10}(1), 22–37.

\phantomsection
\label{csl:82}{Integrating mechanistic and empirical model projections to assess climate impacts on tree species distributions in northwestern North America}. (2017). \textit{Global Change Biology}, \textit{23}(5), 2005–2015. \url{https://doi.org/10.1111/gcb.13570}

\phantomsection
\label{csl:83}{Preferential sampling for presence/absence data and for fusion of presence/absence data with presence-only data}. (2019). \textit{Ecological Monographs}, \textit{89}(3), e01372.

\phantomsection
\label{csl:84}{Integrated models that unite local and regional data reveal larger-scale environmental relationships and improve predictions of species distributions}. (2016). \textit{Landscape Ecology}, \textit{31}(6), 1369–1382.

\phantomsection
\label{csl:85}{Integrated modeling of waterfowl distribution in western Canada using aerial survey and citizen science (eBird) data}. (2021). \textit{Ecosphere}, \textit{12}(10), e03790.

\phantomsection
\label{csl:86}{Identifying and correcting spatial bias in opportunistic citizen science data for wild ungulates in Norway}. (2021). \textit{Ecology and Evolution}, \textit{11}(21), 15191–15204.

\phantomsection
\label{csl:87}{The relative contribution of camera trap technology and citizen science for estimating survival of an endangered African vulture}. (2020). \textit{Biological Conservation}, \textit{246}, 108593.

\phantomsection
\label{csl:88}{Integrating variability in detection probabilities when designing wildlife surveys: a case study of amphibians from south-eastern Australia}. (2012). \textit{Biodiversity and Conservation}, \textit{21}(3), 729–744.

\phantomsection
\label{csl:89}{Estimating species misclassification with occupancy dynamics and encounter rates: A semi-supervised, individual-level approach}. (2022). \textit{Methods in Ecology and Evolution}.

\phantomsection
\label{csl:90}{Accounting for misidentification and heterogeneity in occupancy studies using hidden Markov models}. (2018). \textit{Ecological Modelling}, \textit{387}, 61–69.

\phantomsection
\label{csl:91}{Estimating occupancy from autonomous recording unit data in the presence of misclassifications and detection heterogeneity}. (2022). \textit{Methods in Ecology and Evolution}.

\phantomsection
\label{csl:92}{The current state of citizen science as a tool for ecological research and public engagement}. (2012). \textit{Frontiers in Ecology and the Environment}, \textit{10}(6), 291–297.

\phantomsection
\label{csl:93}{Estimating abundance from repeated presence--absence data or point counts}. (2003). \textit{Ecology}, \textit{84}(3), 777–790.

\phantomsection
\label{csl:94}{Weekend bias in Citizen Science data reporting: implications for phenology studies}. (2013). \textit{International Journal of Biometeorology}, \textit{57}, 715–720.

\phantomsection
\label{csl:95}{Observing the observers: how participants contribute data to iNaturalist and implications for biodiversity science}. (2021). \textit{BioScience}, \textit{71}(11), 1179–1188.

\phantomsection
\label{csl:96}{Rcode and supplementary information for Data Generating Process of Citizen Science Data}. (2022). \textit{Dyrad, Dataset}. \url{https://doi.org/10.5061/dryad.0rxwdbs2m}

\phantomsection
\label{csl:0}Algimantas, P., \& Rasa, U. (2010). {ABUNDANCE AND SEASONAL MIGRATION OF GULLS(LARIDAE) ON THE LITHUANIAN BALTIC SEA COAST}. \textit{Acta Biologica Universitatis Daugavpiliensis}, \textit{10}(2), 147–164.

\phantomsection
\label{csl:0}Augustine, B. C., Royle, J. A., Linden, D. W., \& Fuller, A. K. (2020). {Spatial proximity moderates genotype uncertainty in genetic tagging studies}. \textit{Proceedings of the National Academy of Sciences}, \textit{117}(30), 17903–17912.

\phantomsection
\label{csl:0}Baddeley, A., \& Turner, R. (2005). {Spatstat: an R package for analyzing spatial point patterns}. \textit{Journal of Statistical Software}, \textit{12}(1), 1–42.

\phantomsection
\label{csl:0}Baker, E., Drury, J. P., Judge, J., Roy, D. B., Smith, G. C., \& Stephens, P. A. (2021). {The verification of ecological citizen science data: current approaches and future possibilities}. \textit{Citizen Science: Theory and Practice}, \textit{6}(1).

\phantomsection
\label{csl:0}Bivand, R., Altman, M., Anselin, L., Assun{\c{c}}{\~a}o, R., Berke, O., Bernat, A., \& Blanchet, G. (2015). {Package ‘spdep’}. \textit{The Comprehensive R Archive Network}.

\phantomsection
\label{csl:0}Bonney, R., Cooper, C. B., Dickinson, J., Kelling, S., Phillips, T., Rosenberg, K. V., \& Shirk, J. (2009). {Citizen science: a developing tool for expanding science knowledge and scientific literacy}. \textit{BioScience}, \textit{59}(11), 977–984.

\phantomsection
\label{csl:0}Borowiec, M. L., Dikow, R. B., Frandsen, P. B., McKeeken, A., Valentini, G., \& White, A. E. (2022). {Deep learning as a tool for ecology and evolution}. \textit{Methods in Ecology and Evolution}, \textit{13}(8), 1640–1660.

\phantomsection
\label{csl:0}Brooks, S. P., \& Gelman, A. (1998). {General methods for monitoring convergence of iterative simulations}. \textit{Journal of Computational and Graphical Statistics}, \textit{7}(4), 434–455.

\phantomsection
\label{csl:0}Brym, M. Z., Henry, C., Lukashow-Moore, S. P., Henry, B. J., Van Gestel, N., \& Kendall, R. J. (2020). {Prevalence of monarch (Danaus plexippus) and queen (Danaus gilippus) butterflies in West Texas during the fall of 2018}. \textit{BMC Ecology}, \textit{20}(1), 1–13.

\phantomsection
\label{csl:0}Caradima, B., Schuwirth, N., \& Reichert, P. (2019). {From individual to joint species distribution models: A comparison of model complexity and predictive performance}. \textit{Journal of Biogeography}, \textit{46}(10), 2260–2274.

\phantomsection
\label{csl:0}Chambers, J. (2008). \textit{{Software for data analysis: programming with R}}. Springer Science \& Business Media.

\phantomsection
\label{csl:0}Clare, J. D. J., Townsend, P. A., Anhalt-Depies, C., Locke, C., Stenglein, J. L., Frett, S., Martin, K. J., Singh, A., Van Deelen, T. R., \& Zuckerberg, B. (2019). {Making inference with messy (citizen science) data: when are data accurate enough and how can they be improved?}. \textit{Ecological Applications}, \textit{29}(2), e01849.

\phantomsection
\label{csl:0}Conn, P. B., McClintock, B. T., Cameron, M. F., Johnson, D. S., Moreland, E. E., \& Boveng, P. L. (2013). {Accommodating species identification errors in transect surveys}. \textit{Ecology}, \textit{94}(11), 2607–2618.

\phantomsection
\label{csl:0}Cressie, N. (2015). \textit{{Statistics for spatial data}}. John Wiley \& Sons.

\phantomsection
\label{csl:0}Cruickshank, S. S., B{\"u}hler, C., \& Schmidt, B. R. (2019). {Quantifying data quality in a citizen science monitoring program: False negatives, false positives and occupancy trends}. \textit{Conservation Science and Practice}, \textit{1}(7), e54.

\phantomsection
\label{csl:0}Culbertson, K. A., Garland, M. S., Walton, R. K., Zemaitis, L., \& Pocius, V. M. (2022). {Long-term monitoring indicates shifting fall migration timing in monarch butterflies (Danaus plexippus)}. \textit{Global Change Biology}, \textit{28}(3), 727–738.

\phantomsection
\label{csl:0}Dickinson, J. L., Zuckerberg, B., \& Bonter, D. N. (2010). {Citizen science as an ecological research tool: challenges and benefits}. \textit{Annual Review of Ecology, Evolution, and Systematics}, \textit{41}, 149–172.

\phantomsection
\label{csl:0}Elith, J., \& Leathwick, J. R. (2009). {Species distribution models: ecological explanation and prediction across space and time}. \textit{Annual Review of Ecology, Evolution and Systematics}, \textit{40}(1), 677–697.

\phantomsection
\label{csl:0}Elith*, J., H. Graham*, C., P. Anderson, R., Dud{\'\i}k, M., Ferrier, S., Guisan, A., J. Hijmans, R., Huettmann, F., R. Leathwick, J., Lehmann, A., \& others. (2006). {Novel methods improve prediction of species’ distributions from occurrence data}. \textit{Ecography}, \textit{29}(2), 129–151.

\phantomsection
\label{csl:0}Fahrmeir, L., Kneib, T., Lang, S., \& Marx, B. (2013). {Categorical Regression Models}. In \textit{Regression: Models, Methods and Applications} (pp. 325–347). Springer Berlin Heidelberg. \url{https://doi.org/10.1007/978-3-642-34333-9_6}

\phantomsection
\label{csl:0}Fern{\'a}ndez-i-Mar{\'\i}n, X. (2016). {ggmcmc: Analysis of MCMC samples and Bayesian inference}. \textit{Journal of Statistical Software}, \textit{70}, 1–20.

\phantomsection
\label{csl:0}Gardiner, M. M., Allee, L. L., Brown, P. M. J., Losey, J. E., Roy, H. E., \& Smyth, R. R. (2012). {Lessons from lady beetles: accuracy of monitoring data from US and UK citizen-science programs}. \textit{Frontiers in Ecology and the Environment}, \textit{10}(9), 471–476.

\phantomsection
\label{csl:0}Gelfand, A. E., \& Schliep, E. M. (2018). {Bayesian inference and computing for spatial point patterns}. \textit{NSF-CBMS Regional Conference Series in Probability and Statistics}, \textit{10}, i--125.

\phantomsection
\label{csl:0}Guisan, A., \& Rahbek, C. (2011). {SESAM--a new framework integrating macroecological and species distribution models for predicting spatio-temporal patterns of species assemblages}. In \textit{Journal of Biogeography} (Vol. 38, Number 8, pp. 1433–1444). Wiley Online Library.

\phantomsection
\label{csl:0}G{\'o}mez-Rubio, V., \& Rue, H. (2018). {Markov chain monte carlo with the integrated nested laplace approximation}. \textit{Statistics and Computing}, \textit{28}(5), 1033–1051.

\phantomsection
\label{csl:0}Hijmans, R. J., Van Etten, J., Cheng, J., Mattiuzzi, M., Sumner, M., Greenberg, J. A., Lamigueiro, O. P., Bevan, A., Racine, E. B., Shortridge, A., \& others. (2015). {Package ‘raster’}. \textit{R Package}, \textit{734}.

\phantomsection
\label{csl:0}Hijmans, R. J., van Etten, J., Mattiuzzi, M., Sumner, M., Greenberg, J. A., Lamigueiro, O. P., Bevan, A., Racine, E. B., \& Shortridge, A. (2013). \textit{{Raster package in R}}. Version.

\phantomsection
\label{csl:0}Hoekman, S. T. (2021). {Multi-observer methods for estimating uncertain species identification}. \textit{Ecosphere}, \textit{12}(9), e03648.

\phantomsection
\label{csl:0}Hui, F. K. C., Warton, D. I., \& Foster, S. D. (2015). {Multi-species distribution modeling using penalized mixture of regressions}. \textit{The Annals of Applied Statistics}, \textit{9}(2), 866–882.

\phantomsection
\label{csl:0}Johnston, A., Fink, D., Hochachka, W. M., \& Kelling, S. (2018). {Estimates of observer expertise improve species distributions from citizen science data}. \textit{Methods in Ecology and Evolution}, \textit{9}(1), 88–97.

\phantomsection
\label{csl:0}Johnston, A., Matechou, E., \& Dennis, E. B. (2022). {Outstanding challenges and future directions for biodiversity monitoring using citizen science data}. \textit{Methods in Ecology and Evolution}.

\phantomsection
\label{csl:0}Jones, J. P. G. (2011). {Monitoring species abundance and distribution at the landscape scale}. \textit{Journal of Applied Ecology}, \textit{48}(1), 9–13.

\phantomsection
\label{csl:0}Jongbloed, R. H. (2016). \textit{{Flight height of seabirds: a literature study}}.

\phantomsection
\label{csl:0}Kelling, S., Johnston, A., Bonn, A., Fink, D., Ruiz-Gutierrez, V., Bonney, R., Fernandez, M., Hochachka, W. M., Julliard, R., Kraemer, R., \& others. (2019). {Using semistructured surveys to improve citizen science data for monitoring biodiversity}. \textit{BioScience}, \textit{69}(3), 170–179.

\phantomsection
\label{csl:0}Kelling, S., Johnston, A., Hochachka, W. M., Iliff, M., Fink, D., Gerbracht, J., Lagoze, C., La Sorte, F. A., Moore, T., Wiggins, A., \& others. (2015). {Can observation skills of citizen scientists be estimated using species accumulation curves?}. \textit{PloS One}, \textit{10}(10), e0139600.

\phantomsection
\label{csl:0}Kelling, S., Yu, J., Gerbracht, J., \& Wong, W.-K. (2011). {Emergent filters: automated data verification in a large-scale citizen science project}. \textit{2011 IEEE Seventh International Conference on e-Science Workshops}, 20–27.

\phantomsection
\label{csl:0}Keshavan, A., Yeatman, J. D., \& Rokem, A. (2019). {Combining citizen science and deep learning to amplify expertise in neuroimaging}. \textit{Frontiers in Neuroinformatics}, \textit{13}, 29.

\phantomsection
\label{csl:0}Koch, W., Hogeweg, L., Nilsen, E. B., \& Finstad, A. G. (2022). {Maximizing citizen scientists’ contribution to automated species recognition}. \textit{Scientific Reports}, \textit{12}(1), 1–10.

\phantomsection
\label{csl:0}Lane, M. A., \& Edwards, J. L. (2007). {The global biodiversity information facility (GBIF)}. \textit{Systematics Association Special Volume}, \textit{73}, 1.

\phantomsection
\label{csl:0}Langenk{\"a}mper, D., Simon-Lledo, E., Hosking, B., Jones, D. O. B., \& Nattkemper, T. W. (2019). {On the impact of Citizen Science-derived data quality on deep learning based classification in marine images}. \textit{PloS One}, \textit{14}(6), e0218086.

\phantomsection
\label{csl:0}Lotfian, M., Ingensand, J., \& Brovelli, M. A. (2021). {The Partnership of Citizen Science and Machine Learning: Benefits, Risks, and Future Challenges for Engagement, Data Collection, and Data Quality}. \textit{Sustainability}, \textit{13}(14), 8087.

\phantomsection
\label{csl:0}Madsen, E. K., \& Broekhuis, F. (2020). {Determining multi-species site use outside the protected areas of the Maasai Mara, Kenya, using false positive site-occupancy modelling}. \textit{Oryx}, \textit{54}(3), 395–404.

\phantomsection
\label{csl:0}Marques, O., James, J., \& Barcelos, E. (2018). {Face-It-Up: a scientific app for face processing using mobile devices and machine learning APIs}. \textit{Mobile Multimedia/Image Processing, Security, and Applications 2018}, \textit{10668}, 106680S.

\phantomsection
\label{csl:0}Molinari-Jobin, A., K{\'e}ry, M., Marboutin, E., Molinari, P., Koren, I., Fuxj{\"a}ger, C., Breitenmoser-W{\"u}rsten, C., W{\"o}lfl, S., Fasel, M., Kos, I., \& others. (2012). {Monitoring in the presence of species misidentification: the case of the E urasian lynx in the A lps}. \textit{Animal Conservation}, \textit{15}(3), 266–273.

\phantomsection
\label{csl:0}Njue, N., Kroese, J. S., Gr{\"a}f, J., Jacobs, S. R., Weeser, B., Breuer, L., \& Rufino, M. C. (2019). {Citizen science in hydrological monitoring and ecosystem services management: State of the art and future prospects}. \textit{Science of the Total Environment}, \textit{693}, 133531.

\phantomsection
\label{csl:0}Nugent, J. (2018). {iNaturalist}. \textit{Science Scope}, \textit{41}(7), 12–13.

\phantomsection
\label{csl:0}O'Hara, R. B., \& Sillanp{\"a}{\"a}, M. J. (2009). {A review of Bayesian variable selection methods: what, how and which}. \textit{Bayesian Analysis}, \textit{4}(1), 85–117.

\phantomsection
\label{csl:0}Ovaskainen, O., \& Soininen, J. (2011). {Making more out of sparse data: hierarchical modeling of species communities}. \textit{Ecology}, \textit{92}(2), 289–295.

\phantomsection
\label{csl:0}Perry, J. R., Sumner, S., Thompson, C., \& Hart, A. G. (2021). {‘Citizen identification’: online learning supports highly accurate species identification for insect-focussed citizen science}. \textit{Insect Conservation and Diversity}, \textit{14}(6), 862–867.

\phantomsection
\label{csl:0}Pillay, R., Miller, D. A. W., Hines, J. E., Joshi, A. A., \& Madhusudan, M. D. (2014). {Accounting for false positives improves estimates of occupancy from key informant interviews}. \textit{Diversity and Distributions}, \textit{20}(2), 223–235.

\phantomsection
\label{csl:0}Pollock, L. J., Tingley, R., Morris, W. K., Golding, N., O'Hara, R. B., Parris, K. M., Vesk, P. A., \& McCarthy, M. A. (2014). {Understanding co-occurrence by modelling species simultaneously with a Joint Species Distribution Model (JSDM)}. \textit{Methods in Ecology and Evolution}, \textit{5}(5), 397–406.

\phantomsection
\label{csl:0}Ponti, M., \& Seredko, A. (2022). {Human-machine-learning integration and task allocation in citizen science}. \textit{Humanities and Social Sciences Communications}, \textit{9}(1), 1–15.

\phantomsection
\label{csl:0}Prudic, K. L., Timmermann, B. N., Papaj, D. R., Ritland, D. B., \& Oliver, J. C. (2019). {Mimicry in viceroy butterflies is dependent on abundance of the model queen butterfly}. \textit{Communications Biology}, \textit{2}(1), 1–9.

\phantomsection
\label{csl:0}Renner, I. W., \& Warton, D. I. (2013). {Equivalence of MAXENT and Poisson point process models for species distribution modeling in ecology}. \textit{Biometrics}, \textit{69}(1), 274–281.

\phantomsection
\label{csl:0}Ritland, D. B. (1991). {Revising a classic butterfly mimicry scenario: demonstration of M{\"u}llerian mimicry between Florida viceroys (Limenitis archippus floridensis) and queens (Danaus gilippus berenice)}. \textit{Evolution}, \textit{45}(4), 918–934.

\phantomsection
\label{csl:0}Roberts, D. R., Bahn, V., Ciuti, S., Boyce, M. S., Elith, J., Guillera-Arroita, G., Hauenstein, S., Lahoz-Monfort, J. J., Schr{\"o}der, B., Thuiller, W., \& others. (2017). {Cross-validation strategies for data with temporal, spatial, hierarchical, or phylogenetic structure}. \textit{Ecography}, \textit{40}(8), 913–929.

\phantomsection
\label{csl:0}Rue, H., Lindgren, F., Simpson, D., Martino, S., Teixeira Krainski, E., Bakka, H., Riebler, A., \& Fuglstad, G. A. (2019). {INLA: Full Bayesian analysis of latent Gaussian models using integrated nested Laplace approximations}. \textit{R Package Version}, \textit{19}(03).

\phantomsection
\label{csl:0}Ruiz-Gutierrez, V., Hooten, M. B., \& Grant, E. H. C. (2016). {Uncertainty in biological monitoring: a framework for data collection and analysis to account for multiple sources of sampling bias}. \textit{Methods in Ecology and Evolution}, \textit{7}(8), 900–909.

\phantomsection
\label{csl:0}Santos-Fernandez, E., \& Mengersen, K. (2021). {Understanding the reliability of citizen science observational data using item response models}. \textit{Methods in Ecology and Evolution}, \textit{12}(8), 1533–1548.

\phantomsection
\label{csl:0}Saoud, Z., Fontaine, C., Lo{\"\i}s, G., Julliard, R., \& Rakotoniaina, I. (2020). {Miss-identification detection in citizen science platform for biodiversity monitoring using machine learning}. \textit{Ecological Informatics}, \textit{60}, 101176.

\phantomsection
\label{csl:0}Sauer, J. R., Link, W. A., Fallon, J. E., Pardieck, K. L., \& Ziolkowski, D. J. (2013). {The North American breeding bird survey 1966--2011: summary analysis and species accounts}. \textit{North American Fauna}, \textit{79 (79)}, 1–32.

\phantomsection
\label{csl:0}Sicacha-Parada, J., Pavon-Jordan, D., Steinsland, I., May, R., Stokke, B., \& {\O}ien, I. J. (2021). {A spatial modeling framework for monitoring surveys with different sampling protocols with a case study for bird populations in mid-Scandinavia}. \textit{ArXiv Preprint ArXiv:2104.05751}.

\phantomsection
\label{csl:0}Sicacha-Parada, J., Steinsland, I., Cretois, B., \& Borgelt, J. (2021). {Accounting for spatial varying sampling effort due to accessibility in Citizen Science data: A case study of moose in Norway}. \textit{Spatial Statistics}, \textit{42}, 100446.

\phantomsection
\label{csl:0}Steel, G. D., \& others. (1953). \textit{{Relation between poisson and multinomial distributions}}.

\phantomsection
\label{csl:0}Suzuki-Ohno, Y., Westfechtel, T., Yokoyama, J., Ohno, K., Nakashizuka, T., Kawata, M., \& Okatani, T. (2022). {Deep learning increases the availability of organism photographs taken by citizens in citizen science programs}. \textit{Scientific Reports}, \textit{12}(1), 1–10.

\phantomsection
\label{csl:0}Swanson, A., Kosmala, M., Lintott, C., \& Packer, C. (2016). {A generalized approach for producing, quantifying, and validating citizen science data from wildlife images}. \textit{Conservation Biology}, \textit{30}(3), 520–531.

\phantomsection
\label{csl:0}Tang, B., Clark, J. S., \& Gelfand, A. E. (2021). {Modeling spatially biased citizen science effort through the eBird database}. \textit{Environmental and Ecological Statistics}, \textit{28}(3), 609–630.

\phantomsection
\label{csl:0}Tiago, P., Pereira, H. M., \& Capinha, C. (2017). {Using citizen science data to estimate climatic niches and species distributions}. \textit{Basic and Applied Ecology}, \textit{20}, 75–85.

\phantomsection
\label{csl:0}Tobler, M. W., K{\'e}ry, M., Hui, F. K. C., Guillera-Arroita, G., Knaus, P., \& Sattler, T. (2019). {Joint species distribution models with species correlations and imperfect detection}. \textit{Ecology}, \textit{100}(8), e02754.

\phantomsection
\label{csl:0}Tulloch, A. I. T., Possingham, H. P., Joseph, L. N., Szabo, J., \& Martin, T. G. (2013). {Realising the full potential of citizen science monitoring programs}. \textit{Biological Conservation}, \textit{165}, 128–138.

\phantomsection
\label{csl:0}Ueda, K. (2020). {iNaturalist Research-grade Observations. iNaturalist. org}. \textit{Occurrence Dataset Https://Doi. Org/10.15468/ab3s5x [Accessed via GBIF. Org on 26 May 2020, Https://Www. Gbif. Org/Occurrence/2294559370]}.

\phantomsection
\label{csl:0}Vermeiren, P., Reichert, P., \& Schuwirth, N. (2020). {Integrating uncertain prior knowledge regarding ecological preferences into multi-species distribution models: Effects of model complexity on predictive performance}. \textit{Ecological Modelling}, \textit{420}, 108956.

\phantomsection
\label{csl:0}Vohland, K., Land-Zandstra, A., Ceccaroni, L., Lemmens, R., Perell{\'o}, J., Ponti, M., Samson, R., \& Wagenknecht, K. (2021). \textit{{The science of citizen science}}. Springer Nature.

\phantomsection
\label{csl:0}Wickham, H., Chang, W., \& Wickham, M. H. (2016). {Package ‘ggplot2’}. \textit{Create Elegant Data Visualisations Using the Grammar of Graphics. Version}, \textit{2}(1), 1–189.

\phantomsection
\label{csl:0}Willi, M., Pitman, R. T., Cardoso, A. W., Locke, C., Swanson, A., Boyer, A., Veldthuis, M., \& Fortson, L. (2019). {Identifying animal species in camera trap images using deep learning and citizen science}. \textit{Methods in Ecology and Evolution}, \textit{10}(1), 80–91.

\phantomsection
\label{csl:0}Young, B. E., Dodge, N., Hunt, P. D., Ormes, M., Schlesinger, M. D., \& Shaw, H. Y. (2019). {Using citizen science data to support conservation in environmental regulatory contexts}. \textit{Biological Conservation}, \textit{237}, 57–62.

\phantomsection
\label{csl:0}Yu, H., Jiang, S., \& Land, K. C. (2015). {Multicollinearity in hierarchical linear models}. \textit{Social Science Research}, \textit{53}, 118–136.

\phantomsection
\label{csl:0}Zipkin, E. F., Ries, L., Reeves, R., Regetz, J., \& Oberhauser, K. S. (2012). {Tracking climate impacts on the migratory monarch butterfly}. \textit{Global Change Biology}, \textit{18}(10), 3039–3049.

\phantomsection
\label{csl:0}Zuur, A. F. (2017). {Beginner’s Guide to Spatial, Temporal and Spatial-Temporal Ecological Data Analysis with R-Inla: Using Glm and Glmm Volume I}. \textit{Hightland Statistics Ltd., Sl OCLC}, \textit{973745327}.

\phantomsection
\label{csl:0}de Valpine, P., Paciorek, C., Turek, D., Michaud, N., Anderson-Bergman, C., Obermeyer, F., Wehrhahn Cortes, C., Rodr{\i}guez, A., Temple Lang, D., \& Paganin, S. (2020). {NIMBLE user manual}. \textit{R Package Manual Version 0.9}, \textit{1}.

\phantomsection
\label{csl:0}de Valpine, P., Turek, D., Paciorek, C. J., Anderson-Bergman, C., Lang, D. T., \& Bodik, R. (2017). {Programming with models: writing statistical algorithms for general model structures with NIMBLE}. \textit{Journal of Computational and Graphical Statistics}, \textit{26}(2), 403–413.

\phantomsection
\label{csl:0}{GBIF.Org}. (2022). \textit{{GBIF: The Global Biodiversity Information Facility (year) What is GBIF?}}. The Global Biodiversity Information Facility. \url{https://doi.org/10.15468/DL.H24BP5}

\phantomsection
\label{csl:0}{R Core Team}. (2022). \textit{{R: A Language and Environment for Statistical Computing}}. R Foundation for Statistical Computing. {\url{https://www.R-project.org/}}

\end{document}

