\section{Methods}
\label{sec:Methods}
In this chapter we present the methods that we used to predict the flux of plastic to the ocean (P\textsubscript{O}). Figure \ref{fig:ProcessDiagram} gives a graphical representation of the several steps needed to predict P\textsubscript{O}, which will be explained in the following sections. We first present the data that we collected for the predictor variables (PV) that we used to model the amount of collected waste (MSW\textsubscript{C}) as well as the data including information on the waste collection system and waste treatment (used to predict the MSW\textsubscript{T}) in section \ref{sec:DataSources}. We further present the correlations and trends among the PVs and with the target variable in section \ref{sec:TrendsAndCorrelations}. This is followed by a description of the modification of said data (to allow use in the machine learning models (MLM)) in section \ref{sec:DataModification}. Then, we describe the models that we used to predict the MSW\textsubscript{C}, including the model optimization (section \ref{sec:MLMOptimization}), the assessment of the PV importance (section \ref{sec:PredictorVariableSelection}) and the selection of the best model (section \ref{sec:SelectionOfBestModel}). Following, we introduce our waste flow model (WFM), including the calculation of the total amount of MSW (MSW\textsubscript{T}, section \ref{sec:CalculatingMSWT}), plastic waste that enters the environment (P\textsubscript{E}, section  \ref{sec:CalculatingPlasticEnteringEnvironment}) and plastic waste that enters the ocean (P\textsubscript{O}, section \ref{sec:ProbabilityOfTransport}).
\begin{figure*}[htbp]
	\renewcommand{\arraystretch}{1.2}
	\centering
	\includegraphics[width=\textwidth]{img/process_diagram.pdf}
	\caption{Summary of the processes used to estimate P\textsubscript{O}. Each box represents one major step and indicates which section describes this process. The arrows indicate which data or models were used in the process. The column \textit{Data preparation} describes how the data was obtained and modified to be used in the MLM. The column \textit{machine learning models} describes how the MLM were optimized, how they were used to select the optimal PVs and how the the MSW\textsubscript{C} was estimated. Finally, the column \textit{waste flow model} describes how we estimated MSW\textsubscript{T} and P\textsubscript{T} using information on the waste collection and composition data. The WFM used to estimate P\textsubscript{E} and P\textsubscript{O} is further explained in Figure \ref{fig:WFM}.}
	\label{fig:ProcessDiagram}
\end{figure*}

\subsection{Description of the Data}
\label{sec:DataSources}
Here we present all the data we obtained to predict the plastic input to the ocean (compare Fig. \ref{fig:ProcessDiagram}). We begin with the description of the available data for MSW\textsubscript{C}, followed by the description of the PV data used in the MLM to predict the MSW\textsubscript{C}. Then we present the data of the waste collection system, the waste composition and the waste disposal, which was used to predict the total amount of plastic waste (P\textsubscript{T}) from the MSW\textsubscript{C}. This data is summarized in table \ref{tab:PredictorVariables}. Finally, we present the geospatial data that we use to predict P\textsubscript{E} on a global raster as well as in the WFM and for the ocean current advection model. A detailed overview of the data presented in this section can be found in table \ref{tab:DataCompleteness} in the appendix. This table indicates for each variable how many data points were available per country.\\
	\subsubsection{MSW\textsubscript{C} data}
	We obtained data for our target variable MSW\textsubscript{C} from the UN database. The dataset includes 1'544 data points for 125 countries in the time span 1990-2015 (compare table \ref{tab:PredictorVariables}). No country has a complete dataset for the whole time period (1990-2015, 26 years). However, for four countries (Austria, Luxembourg, Poland and Switzerland) only the data for 2015 is missing. A total of 27 countries has data available for at least 20 years and 64 countries provide data for more than 10 years. 23 countries provide data for from two up to five years and six countries (Angola, Bhutan, Costa Rica, Peru, Sri Lanka and Uganda) provide data for only one year. The data availability is highest in the years 2000-2012 where more than 65 countries have data available.

	\subsubsection{PV data}
	We obtained the data for the PVs from the UN database and the World Bank database. The data contains information for 195 countries (officially recognized by the UN) as well as for 36 dependencies and autonomous regions (e.g. Hong Kong, Falkland Islands). For simplicity, these countries, dependencies, and autonomous regions will further on simply be referred to as countries. Table \ref{tab:PredictorVariables} gives a summary of each dataset. Data for the PVs \textit{total population} and \textit{urban population} is available since 1950 and includes projections until 2050. These two datasets cover all 231 countries. The dataset \textit{young population} covers fewer countries (221) and over a shorter time scale from 1960 until 2017. \textit{GDP} covers the time period 1960-2015 and includes data for 202 countries. Data on \textit{greenhouse gas} and \textit{electricity access} was available from 1970-2012 and for 195 ans 203 countries respectively. Data for F\textsubscript{recycled} is available for 102 countries for the years 1990 and 1995-2015.
		
	\subsubsection{Data of waste collection, composition and disposal}
	We obtained the data on the fraction of population with access to a waste collection system (F\textsubscript{collection}), the fraction of plastic in MSW\textsubscript{C} and the fractions of MSW\textsubscript{C} that is deposited in landfills (F\textsubscript{landfills}) or dumps (F\textsubscript{dumps}) from the UN database and from Hoornweg\&Bhada-Tata \cite{Hoornweg2012}. Data for F\textsubscript{landfill} is available for 102 countries in the time period 1990-2015. Data for F\textsubscript{recycled} is only available for 41 countries and in the years 1990 and 1995-2013. Data for F\textsubscript{landfills} and F\textsubscript{dumps} is only available for 74 countries in the year 2012. Data for the GNI is available from the World Bank for 204 countries in the time period 1960-2016.
	
	\begin{table*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\rowcolors{37}{color2!10}{}
		\caption{Description and sources of the target variable MSW\textsubscript{C} of the MLM, all PVs used in the MLM to predict MSW\textsubscript{C} for each country and year and the variables used to estimate MSW\textsubscript{T}. \textit{N countries} indicates for how many countries at least one data point was available. \textit{N data points} indicates how many data points were available (including the represented fraction of the maximum possible data points for this time period and number of countries).}
		\begin{tabularx}{\textwidth}{p{2.9cm} >{\centering\arraybackslash}p{1cm} p{1.6cm} p{1.5cm} >{\centering\arraybackslash}p{2.3cm} X}
		\toprule
		\textbf{Data}&\textbf{Source}&\textbf{Time scale}&\textbf{N countries}&\textbf{N data points}&\textbf{Description}\\ \toprule
		\multicolumn{6}{l}{\textbf{Target variable}}\\ \midrule
		MSW\textsubscript{C} &\cite{UnitedNations2018MSWcollected}	&1990-2015	&125	& 1'544 (48\%) &Amount of MSW that was collected in a year.\\ \midrule
		\multicolumn{6}{l}{\textbf{PVs used in the MLM to predict MSW\textsubscript{C}}}\\ \midrule
		Total population	&\cite{UnitedNations2018Pop}	&1950-2050	&231	& 23'331 (100\%) &Number of people living in the country.\\
		Urban population	&\cite{UnitedNations2018Urb}	&1950-2050	&231	& 23'331 (100\%) &Number of people living in urban areas.\\
		Young population 	&\cite{TheWorldBank2014Young}	&1960-2017	& 195	& 11'241 (99\%) &Fraction of the total population aged 0-14 years.\\
		Electricity access	&\cite{TheWorldBank2018Ele}	&1990-2014	&204	& 4'931 (83\%) &Fraction of households with access to electricity.\\
		GDP					&\cite{TheWorldBank2018GDP}	&1960-2016	&202	& 8'913 (77\%) &Gross domestic product per country in USD, current conversion rate.\\
		Energy consumption 	&\cite{TheWorldBank2018Ene}	&1990-2015	&173	& 3'644 (78\%) &Energy consumption per person (in kg oil equivalents).\\
		Greenhouse gas		&\cite{TheWorldBank2018GHG}	&1970-2012	&199	& 8'060 (94\%) &Total emissions of greenhouse gases in CO\textsubscript{2} equivalents.\\
		F\textsubscript{recycled}	&\cite{UnitedNations2018MSWrecycled}	&1990, \newline 1995-2013	& 78	& 916 (59\%) &Fraction of MSW\textsubscript{C} that is recycled.\\
		\midrule
		\multicolumn{6}{l}{\textbf{Variables used in to estimate MSW\textsubscript{T}}}\\ \midrule
		GNI		&\cite{TheWorldBank2018GNI}	&1960-2016	&204	& 8'901 (70\%) &Gross national income per country in USD, current conversion rate.\\
		F\textsubscript{collection} &\cite{UnitedNations2018MSWcollection}	&1990-2015	&102	& 958 (36\%) &Fraction of MSW\textsubscript{T} that is collected at the households. \\
		F\textsubscript{plastic} &\cite{Hoornweg2012}		&2012		&103	& 103 &Fraction of plastic in MSW\textsubscript{C}.\\
		F\textsubscript{landfills}, F\textsubscript{dumps}	&\cite{Hoornweg2012}		&2012		&74		& 74 &Fraction of MSW\textsubscript{C} that ends in either landfills or dumps.\\ \bottomrule
		\end{tabularx}
		\label{tab:PredictorVariables}
	\end{table*}
	
	\subsubsection{Geospatial data}
	\textbf{Population density}\\
	We obtained the data on the population density from the Socioeconomic Data and Applications Center (SEDAC). It is available as an ArcInfo GRID file for the years 1990 and 1995 (in a 2.5 arc-minutes resolution) or as a GeoTiff-image for the years 2000, 2005, 2010 and 2015 (in a 30 arc-seconds resolution) (table \ref{tab:GeospatialData}). We processed the data to create data of the same grid size and with continuous information for the years 1990-2015. To this end, the data was loaded into RStudio using the \textit{raster} function from the \textit{rgdal} package \cite{Roger2018} and then converted to a 0.1x0.1° raster (6 arc-minutes) using the \textit{aggregate} function from the \textit{raster} package \cite{Hijmans2017}. The value of each aggregated cell was calculated using the mean of all original cells. Some aggregated cells (at the shore line) of the 2.5 arc-minute rasters were therefore outside the scope of the target raster (e.g. ocean cells). The aggregated cells were compared to the rasterized coastline and cells that laid outside the target raster (e.g. in the ocean) were dropped. To create continuous data, the population was then interpolated for the missing years. This was done using the growth rate of each country's population as obtained from the predictor variable \textit{total population}. However, the urbanization in this time period was neglected and the population thus increases linearly in every cell.\\
	
	\textbf{River shape}\\
	We obtained a shape file including 342 river centrelines from NaturalEarth (table \ref{tab:GeospatialData}). This data is used in the WFM and includes all major rivers at the global scale. Lakes are represented with their centrelines and not with their full area. Intermittent rivers (large, seasonal rivers) are not included in this dataset. We created a raster of all rivers using the \textit{rasterize} function of the \textit{raster} package \cite{Hijmans2017} for further use in the WFM. \\
	
	\textbf{Ocean currents}\\
	We obtained the ocean currents data from the GlobCurrent project \cite{GlobCurrentConsortium2015}. The dataset includes daily averaged satellite-based estimates of the surface Eulerian total current from 1993 to 2017 on a 0.25x0.25° global grid (table \ref{tab:GeospatialData}). The Eulerian total current, representing the sum of geostrophic and Ekman current components, was chosen over the simpler pure geostrophic or Ekman current datasets as it was found to represent drifting particle distributions more realistic \cite{GlobCurrentConsortium2015}. However, a known issue of this dataset is that current velocities in local areas close to the coast can be erroneous \cite{GlobCurrentConsortium2015}. Furthermore, the Mediterranean is not included in the dataset and particles released there will thus not be advected. A graphical representation of the data can be seen in Figure \ref{fig:Current_data}.
	\begin{figure*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\centering
		\includegraphics[]{img/current_data.pdf}
		\caption{Visualization of the GlobCurrent surface Eulerian total current data in the South China Sea and the western North Pacific as used in the PARCELS simulations, including the current velocity as the background layer and indicating current streamlines in green. Visualization created using the Syntool visualization website and manually edited \cite{OceanDataLab}.}
		\label{fig:Current_data}
	\end{figure*}
	\begin{table*}[htbp]
	\renewcommand{\arraystretch}{1.2}
	\rowcolors{2}{}{color2!10}
	\caption{Description and sources of the geospatial data. The population density is used in the prediction of MSW\textsubscript{C} on a 0.1x0.1° raster. The river dataset is used in the calculation of plastic transport to the oceans. The ocean current data is used to simulate the advection of plastic particles in the ocean.}
	\begin{tabularx}{\textwidth}{p{2.9cm} >{\centering\arraybackslash}p{1cm} p{1.6cm} X X}
		\toprule
		\textbf{Data}&\textbf{Source}&\textbf{Time scale}&\textbf{Type}&\textbf{Description}\\ \toprule
		Population density &\cite{SEDAC2005}	&1990, 1995	&Raster with a resolution of 2.5 arc-minutes	&Population density.\\
		Population density &\cite{SEDAC2017}	&2000, 2005, 2010, 2015	&Raster with a resolution of 30 arc-seconds	&Population density.\\
		Rivers 		&\cite{NaturalEarth}	& & Shape file of lake and river centres with a resolution of 1:10\textsuperscript{6} & Polyline of major river and lake centres.\\
		Ocean currents & \cite{GlobCurrentConsortium2017} & 1993-2017 & NetCDF file with a resolution of 0.25x0.25°. & Total Eulerian currents (sum of geostrophic and Ekman components). \\ \bottomrule
	\end{tabularx}
	\label{tab:GeospatialData}
\end{table*}

In this section we described the availability and completeness of the individual datasets. In the following section we present the correlations among the datasets.\\
\vfill\null

\subsection[Trends and Correlations among the Predictor \newline Variables and with the Target Variable]{Trends and Correlations among the Predictor Variables and with the Target Variable}
\label{sec:TrendsAndCorrelations}
%Show the plots of var v var and var v MSW. (lin_fit and 95percentile.pdf) \\
\begin{figure*}[htbp]
	\renewcommand{\arraystretch}{1.2}
	\centering
	\includegraphics[width=\textwidth]{img/correlation.pdf}
	\caption{Pearson and Spearman correlations between MSW\textsubscript{C} and the predictor variables as well as among the predictor variables. We expect that highly correlated variables are removed in the variable selection procedure explained in section \ref{sec:PredictorVariableSelection}.}
	\label{fig:Correlation}
\end{figure*}
In this section we present trends and correlations in the datasets that were described in the previous section. We first describe the correlations among the PVs and our target variable and compare this to our expectations as described in the introduction. Then we show correlation among the PVs themselves.\\

\textbf{Correlations among MSW\textsubscript{C} and the PVs}\\
Here we present the Pearson and Spearman correlations of MSW\textsubscript{C} with the PVs (Fig. \ref{fig:Correlation}). The Pearson correlation measures the relationship between linearly related variables while the Spearman correlation can be used to assess monotonic or exponential relationships.\\
The Pearson correlation (Fig. \ref{fig:Correlation} A) of the PVs \textit{energy consumption} (0.52) and  \textit{GDP} (0.52) are \textgreater0.5. \textit{Young population} (-0.31) is the only negative correlation. \textit{F\textsubscript{recycled}} (0.24) has the correlation closest to zero. The Spearman correlation (Fig. \ref{fig:Correlation} B) shows similar values as the Person correlation. Here, \textit{electricity access} is the only PV whose difference to the Pearson correlation is \textgreater0.07. This can be explained by looking at the scatterplot of \textit{electricity access} and MSW\textsubscript{C} in Figure \ref{fig:Scatterplots} in the appendix. It can be seen that \textit{electricity access} the data is not linearly distributed and shows an exponential trend which is captured by the Spearman correlation.\\

\textbf{Correlations among the PVs}\\
The Pearson correlation (Fig. \ref{fig:Correlation} A) reveal that the correlations among the PVs cover a wider range (0.18 - 0.74) than with the target variable MSW\textsubscript{C}. Correlations \textgreater0.5 can be found between \textit{young population} and \textit{electricity access} (-0.74), between \textit{greenhouse gas} and \textit{energy consumption} (0.61), between \textit{energy consumption} and \textit{urban population} (0.57) and between \textit{urban population} and \textit{electricity access} (0.55). We expect that if two variables are highly correlated, one will be removed in the model optimization process (described in section \ref{sec:PredictorVariableSelection}). Highly correlated variables include similar relationships towards MSW\textsubscript{C} thus potentially redundant information.\\
The Spearman correlation (Fig. \ref{fig:Correlation} B) among the PV cover a slightly smaller range (0.17 - 0.63 ) than the Pearson correlations (0-18 - 0.74) and a slightly larger range than with MSW\textsubscript{C} (0.27 - 0.58). Compared to the Pearson correlations, \textit{electricity access} shows a stronger Spearman correlation for all PVs except \textit{urban population} and \textit{young population}. When looking at the scatterplots of these PVs \ref{fig:ScatterplotsPV}, a linear relationship can be seen, which explains the higher Pearson correlation. Furthermore, the Spearman correlations between \textit{GDP} and \textit{greenhouse gas} (0.43) and between \textit{urban population} and \textit{young population} (0.17) are lower respectively higher than the Pearson correlations (0.29 and -0.44 respectively). The difference in the Spearman and Pearson correlations among the other PVs is \textless0.15.\\
\vfill\null

To further visualize the relationships between the PVs, we created a dendrogram that uses the correlation among the PVs as  measure of similarity (Fig. \ref{fig:CorrelationDendrogram} in the appendix). If can be seen that the PVs \textit{electricity access} and \textit{young population} as well as \textit{energy consumption} and \textit{greenhouse gas} are grouped in both the Pearson and Spearman dendrogram. Additionally, \textit{GDP} and MSW\textsubscript{C} are also grouped. Again, we expect that, if two PVs are grouped, one will be removed in the model optimization process.\\
In the next section we describe how we prepared the data for the use in the MLM.

\subsection{Data Preparation for the MLM}
\label{sec:DataModification}
The obtained data as described above is often incomplete regarding the number of countries or years that are represented as well as the number of data points (compare table \ref{tab:PredictorVariables}). As the MLM and the WFM need complete PV datasets, the data had to be prepared for the modelling steps (compare Fig. \ref{fig:ProcessDiagram}). We first reduced and simplified redundant data and then inter-and extrapolated missing data. This increased the data completeness by filling the gaps. Following, we rasterized the data, which is needed for the geospatial analyses in the WFM.
%Only the datasets including population data cover a wide range of \textgreater30 years with a high data completeness (\textgreater99\%). All the other datasets presented in table \ref{tab:PredictorVariables} have more missing values and cover shorter time periods. The lower data completeness can be explained with a poor data availability in the years \textless1990. Also our target variable MSW\textsubscript{C} is only available from 1990-2015.
	\subsubsection{Reduction and simplification of redundant data}
	Here we describe how we reduced and simplified the datasets to increase the data completeness. Following measures were performed:
	\begin{itemize}[noitemsep,topsep=0pt]
		\item We grouped small islands states in Oceania (e.g. Federated States of Micronesia, Marshall Islands, N=16) as well as in the Caribbean (e.g. Jamaica, Cayman Islands, N=17) and applied the mean values (compare table \ref{tab:DroppedCountries}). We thus assume that they are similar in their MSW production rates.
		\item We dropped countries that did not provide GDP data (e.g. Niger, North Korea). This includes 28 countries as presented in table \ref{tab:DroppedCountries}.
		\item We cropped all datasets to the years 1990-2015, as these years contain substantially more data than in the period 1960-1989. Data for 2016 and 2017 was often not yet published.
	\end{itemize}
	These measures reduced the number of countries that we used to model the MSW\textsubscript{C} from 231 to 173. The dropped countries represent 0.78\% of the world's population. Table \ref{tab:DroppedCountries} gives an overview of the dropped and grouped countries. Together with the reduction of the time period, we reduced the number of needed data points from 15'708 (231 countries and 68 years from 1950-2017) to 4'498 (173 countries and 26 years from 1990-2015). This increased the data completeness of the modified datasets (compare table \ref{tab:DataModification}). As the MLM require complete datasets for the PVs, and these measures did not achieve 100\% data completeness, we had to estimate the missing data.\\
	We found two cases of missing data. First, countries that provided PV values for several years but the data does not cover the whole time period. This missing data was inter or extrapolated as described in the next subsection. Secondly, countries that did not provide data for a PV. In this case, the data was estimated based on similar countries. This process is described section \ref{sec:Clustering}.
	
	\begin{table*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\rowcolors{2}{}{color2!10}
		\centering
			\caption{Countries, dependencies and autonomous regions that were dropped or grouped from the analysis. Dropped countries did not provide GDP data which was needed for creating clusters that allow the estimation of missing data (section \ref{sec:MissingData}). We grouped the Caribbean and Oceania islands because we assume that they behave similarly regarding the MSW\textsubscript{C} production.}
			\begin{tabularx}{\textwidth}{X X X X X X}
			\toprule
			\multicolumn{6}{l}{\textbf{Dropped due to missing data}}\\ \toprule
			Anguilla & Caribbean Netherlands  & Cook Islands & Curaçao & Falkland Islands & French Guiana\\
			Gibraltar & Guadeloupe & Holy See & Kyrgyzstan & Martinique & Mayotte\\
			Montserrat & Niger & Niue & North Korea & Palestine & Réunion \\
			Saint Helena & Saint Martin & Saint Pierre and Miquelon & Sint Maarten & Tokelau & Turks and Caicos Islands \\
			Tuvalu & United states Virgin Islands & Wallis and Futuna Islands & Western Sahara & & \\
			\toprule
			\multicolumn{6}{l}{\textbf{Grouped to "Caribbean Islands"}}\\ \hline 
			Antigua and Barbuda & Aruba & Bahamas & Barbados & British Virgin Islands & Cayman Islands \\
			Cuba & Dominica & Dominican Republic & Grenada & Haiti & Jamaica \\
			Puerto Rico & Saint Kitts and Nevis & Saint Lucia & Saint Vincent and the Grenadines & Trinidad and Tobago & \\
			\toprule
			\multicolumn{6}{l}{\textbf{Grouped to "Oceania Islands"}}\\ \hline 
			American Samoa & Federated States of Micronesia & Fiji & French Polynesia & Guam & Kiribati \\
			Marshall Islands  & Nauru & New Caledonia & Northern Mariana Islands & Palau & Papua New Guinea \\
			Samoa & Solomon Islands & Tonga &	Vanuatu &  & \\ \bottomrule
		\end{tabularx} 
		\label{tab:DroppedCountries}
	\end{table*}

	\subsubsection{Inter- and extrapolation of partially missing data}
	\label{sec:MissingData}
	To estimate the missing data of countries where more than 2 data points were available in the years 1990-2015, we calculated the missing values using either a linear or polynomial fit. A linear fit was used for the datasets \textit{electricity access}, \textit{energy consumption}, \textit{greenhouse gas}, F\textsubscript{collection} and \textit{F\textsubscript{recycled}} as these PVs show a linear trend in the available data (compare table \ref{tab:DataModification}). We used a 2nd order polynomial fit for the GDP as the available data suggests a non-linear trend (compare table \ref{tab:DataModification}). We did not allow the calculated data from the linear fit of \textit{F\textsubscript{collection}} to fall below the respective continental minimum as presented by Hoornweg\&Bhada-Tata \cite{Hoornweg2012} (applied values for Africa: 18\%, Asia: 50\% Europe: 60\% and South America: 11\%). Missing data in the the datasets F\textsubscript{dumps} and F\textsubscript{landfills} was replaced by the mean values according to the countries income class (high-income: \textgreater 12'476 USD/year, upper middle-income: 4'046-12'475 USD/year, lower middle-income: 1'026-4'035 USD/year and low-income: \textless 1'025 USD/year) as presented in Hoornweg\&Bhada-Tata \cite{Hoornweg2012}.

	\subsubsection{Estimation for countries with completely missing data}
	\label{sec:Clustering}
	To get an estimate of each variable for countries where no data was available to extrapolate (compare tables \ref{tab:DataModification} and \ref{tab:DataCompleteness}), we replaced the missing data with the mean from countries with a similar GDP. Therefore, we sorted all countries into clusters according to their similarity in their GDP trajectory (Fig. \ref{fig:Clusters}). Countries where no data was available to run a linear or exponential fit were given the mean values of the corresponding cluster (compare table \ref{tab:DataModification}). This was performed on the datasets \textit{electricity access}, \textit{energy consumption}, \textit{greenhouse gas}, \textit{F\textsubscript{collection}}, \textit{F\textsubscript{plastic}} and \textit{F\textsubscript{recycled}}.
	%Even though GDP has been found to be decoupled from MSW production \cite{Hoornweg2015,Eriksson2005,Ausubel2008}, it makes sense to group these countries and focus on the other predictor variables, as the countries, where decoupling is observed, are often highly industrialized and are thus expected to behave similarly in terms of waste production. Additionally, GDP was the only variable where enough data was available to run the cluster algorithm. Clustering based on population size or geographic location might result in groups of very inhomogeneous countries.\\
	The clustering was done using the \textit{kmlShape} package \cite{Christophe2016} version 0.9.5. This package calculates K-means clusters for longitudinal data using shape-respecting distance and shape-respecting means. The ideal number of clusters was determined by minimizing the sum of differences between the trajectory of each country and the respective cluster trajectory mean, running a maximum of 15 clusters (compare Fig. \ref{fig:ClusterError} in the appendix). The optimum result was obtained by creating 7 clusters and is presented in Figure \ref{fig:Clusters}. It can be seen that the clusters 1-6 include countries from at least two continents and only the cluster 7 include solely the European countries Liechtenstein and Monaco (compare also table \ref{tab:DataCompleteness}). Table \ref{tab:DataModification} indicates the number of counties for which the data was estimated using this approach.\\
	The now complete datasets are ready to be used in the MLM, but the geospatial analyses in the WFM require that we rasterize the data. This is described in the next section.	
	\begin{figure*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\centering
		\includegraphics[width=\textwidth]{img/clusters.pdf}
		\caption{Countries clustered by their GDP as described in section \ref{sec:Clustering}. The clusters are used to estimate data of countries that did not provide any data and did thus not allow inter-and extrapolation. The numbers in brackets indicate the average GDP per cluster. Grey countries are not included in the calculations due to missing GDP values (compare the tables \ref{tab:DroppedCountries} and \ref{tab:DataCompleteness}). Note that cluster 7 includes only the countries Monaco and Lichtenstein and is thus not visible on the global map. Table \ref{tab:DataCompleteness} gives detailed information on the GDP and cluster of each country.}
		\label{fig:Clusters}
	\end{figure*}

	\subsubsection{Rasterization of Country Level Data}
	\label{sec:Rasterization}
	After performing the steps described above, the datasets presented in table \ref{tab:PredictorVariables} now have a data completeness of 100\% and are ready to be used in the MLM. For the WFM, the data needs to be rasterized to allow geospatial analyses. Here we describe how we created a 0.1x0.1° raster for each dataset as presented in table \ref{tab:PredictorVariables} for each year in the period 1990-2015.\\
	We obtained a spatial polygon file of the administrative boundaries for each country using the \textit{getData} function from the \textit{raster} package \cite{Hijmans2017}. The data of each dataset was assigned to each cell of the 0.1x0.1° raster whose cell-centre lays in the corresponding spatial polygon. Due to its small area of just 2 km\textsuperscript{2}, the country of Monaco is not represented in any raster as no cell has its centre inside Monaco's boundaries. Thus, the data of Monaco was added manually to the cell containing the city centre of Monaco in each raster. Therefore, the country France "lost" a cell to Monaco. This process was however not needed for other city states or micro states such as Singapore or Kuwait.\\
	%\subsubsection{Linking each land cell to an ocean cell}
	%As the plastic input to the ocean is locally disconnected from the place of waste production, each land and river cell had to be linked to the closest river or ocean cell. This was done by letting a focal window run along the shoreline and determining the nearest shore cell for each land cell based on the \textit{distance} function from the \textit{raster} package \cite{Hijmans2017}. The same was done for the river cells, which were in turn linked to the ocean cell at their river mouth. Thus, each land cell in 100 km from the shore or 20 km from a river was linked to an ocean cell. The plastic input into each ocean cell could now be calculated by summarizing the values of each of its linked land cells.
		
With this, the data is now prepared to be used in the MLM and WFM. In the next section we describe how we used the MLM to predict MSW\textsubscript{C} using a selection of PVs and section \ref{sec:WFM} describes the WFM.

\subsection{Modelling MSW\textsubscript{C}}
\label{sec:ModellingMSWProduction}
After describing the acquisition and preparation of the PVs we now present the modelling approach used to predict MSW\textsubscript{C}. As the observed trends between the PVs and MSW\textsubscript{C} are not always linear (compare section \ref{sec:TrendsAndCorrelations}) we chose to use both NNET and RF models to predict MSW\textsubscript{C}. To avoid overfitting, we run tests to find the best compromise between the least complex architecture, a minimal amount of PVs and a minimal error in the predictions. Less complex models do not only prevent overfitting but were also found to be better in predicting unseen data \cite{Ozesmi2006}. Furthermore we asses the stability of the models to the size of the training data and biases in the training data. We use two different models, including a neural network (NNET) and a random forest (RF) which allows a more robust estimate as the model results can be averaged. Both models use the available MSW\textsubscript{C} values as target variable and the PVs as described in section \ref{sec:DataModification} and table \ref{tab:PredictorVariables}. We first describe how we optimized the architecture of the models. Following, we describe how we selected the best PVs and then show how we assessed the model stability. In the end of this section we describe how we selected the best NNET and RF model that will then be used to predict MSW\textsubscript{C}.

	\subsubsection{Model architecture optimization}
	\label{sec:MLMOptimization}
	We used a individual procedure for the NNET and RF model to optimize the architectures. In both cases we used all points where MSW\textsubscript{C} data was available (N=1'394). We selected the corresponding values of the PVs and normalized them. The target variable MSW\textsubscript{C} was also normalized. This was done to prevent a bias in the models \cite{Ozesmi2006}. We then performed following steps:\\
	
	\textbf{Architecture optimization of the NNET}\\
	The NNET was created using the R package \textit{neuralnet} \cite{Fritsch2016} version 1.33. We tested several network architectures by training ten models on the same, randomly chosen 90\% of the scaled data. To evaluate the model we performed a 10-fold cross validation (CV) and calculated the mean squared error (MSE). This process was repeated ten times on newly randomized data samples to avoid a bias in the training data. Thus, a total of 100 individual models were created. We retained the model with the lowest MSE and compared it to the best model of the other model architectures (Fig. \ref{fig:NNETArchitecture}). Our optimal architecture for the NNET model consists of one layer with five nodes (MSE:0.03). While adding more nodes to the layer does not decrease the MSE, removing nodes results in a loss of accuracy. While some architectures consisting of 2 layers (e.g. three and two nodes per layer, Fig. \ref{fig:NNETArchitecture_all}) were equally successful in prediction, we preferred the one layer model as to limit the model complexity. A summary of all layer and node combinations can be found in Figure \ref{fig:NNETArchitecture_all} in the appendix.\\
	
	\textbf{Architecture optimization of the RF}\\
	The RF was created using the R package \textit{randomForest} \cite{Breiman2018} version 4.6-14. The RF creates a "forest" of 500 unique decision trees and averages the outcome to obtain the best prediction. We chose 500 trees as increasing this number did not decrease the MSE (Fig. \ref{fig:RF_error_ntree}) \cite{Svetnik2003}. Each tree has a random element in its creation and only selects a subset of the PVs in each branch. We optimized the size of this subset by using the \textit{tuneRF} function, which calculates the out of bag error for each try \cite{Breiman2018}. We run the \textit{tuneRF} function 100 times using different randomly sampled 90\% of the data to avoid a bias in the training data \cite{Svetnik2003}. For our data, a \textit{mtry} value of two results in the smallest out of bag error in 96\% of the tries while in 4\% using three variables is preferred. To reduce complexity and prevent overfitting \cite{Svetnik2003}, we thus set \textit{mtry} to two.
	\begin{figure}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\centering
		\includegraphics[width=\columnwidth]{img/NNET_architecture.pdf}
		\caption{Summary of the average MSE for the different one layer model architectures for the NNET. The first number indicates the amount of PVs, the second number shows the number of nodes in the hidden layer and the third number stands for the number of outputs. To avoid overfitting of the NNET we chose the simplest architecture that has a low MSE.}
		\label{fig:NNETArchitecture}
	\end{figure}
	
	\subsubsection{Predictor variable selection}
	\label{sec:PredictorVariableSelection}
	Besides using a simple architecture, using a minimal number of PV can help to prevent overfitting \cite{Ozesmi2006, Svetnik2003}. We ran several tests to select which PVs are most important in the model creation. We use the following techniques to assess the importance of each PV: removal, random permutation and univariate runs. Following we describe these procedures:\\
	\textbf{Removal:} We run both MLM while sequentially removing one parameter at a time. By removing one parameter from the training data we assess the importance of this variable in the model creation \cite{Ozesmi2006}. For each removed PV, we trained ten models on the same randomly sampled 90\% of the data. This was repeated nine more times using newly sampled training data to avoid a bias in the training data. Thus, a total of 100 models were created. We then remove the parameter that has the least effect on the model performance (MSE). This process was repeated with the remaining PVs until only one PV was left.\\
	\textbf{Random permutation:} Analogue to the removal process we run both models while randomly permuting one parameter at a time. This means that for one PV, the values were randomly replaced with values from within the same PV. If the PV is important in the model creation this is expected to create larger errors than for less important PVs \cite{Svetnik2003}. For each randomly permuted PV, we trained ten models on the same randomly sampled 90\% of the data. This was repeated nine more times using newly sampled training data to avoid a bias in the training data. We then remove the parameter that has the least effect on the model performance (MSE). This process was repeated with the remaining PVs until only one PV was left.\\
	\textbf{Univariate:} We ran the model using only one PV at a time and calculated the MSE. Here, we expect that PVs that are important for the model creation have lower MSE values \cite{Ozesmi2006}.\\
	
	We then compared the effect of the removal, random permutation and univariate runs on the MSE. Figure \ref{fig:VariableSelection} shows the example of PV removal in the RF model. As soon as removing a PV created a big increase in the MSE, we assumed that this PV, together with the remaining PVs, is crucial in the model creation \cite{Ozesmi2006}. In this case, removing the PVs, \textit{F\textsubscript{recycled}}, \textit{GDP}, \textit{energy consumption} or \textit{electricity access} has only a small effect on the MSE. If \textit{young population} is removed however, the effect on the MSE is stronger. Thus, we selected the three PVs \textit{young population}, \textit{greenhouse gas} and \textit{urban population} as essential for the model.\\
	A summary of the PV selection (including both models and all selection processes as described above) can be seen in table \ref{tab:VariableSelection}. While no PV was always the most or least important, clear trends can be seen. The PVs \textit{urban population} and \textit{greenhouse gas} are the most important in both the average rank (1.2 and 2.2 respectively) as well as the median rank (1 and 2 respectively). These two PVs were always among the three most important PVs throughout all selection types. The PVs \textit{young population} and \textit{energy consumption} have an average rank of 4.0, but the median of latter is higher (3.5 and 5 respectively). The parameter removal and parameter shuffling processes both favour \textit{young population} and only the univariate selection process favours \textit{energy consumption}. Thus, \textit{young population} was preferred. F\textsubscript{recycled} was usually unimportant and thus has the lowest average and median rank (6.5 in both). As a result of these tests we selected \textit{young population}, \textit{greenhouse gas} and \textit{urban population} as our most important PVs. They were used in the prediction of MSW\textsubscript{C} while the others were dropped.
	
	\begin{figure*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\centering
		\includegraphics[width=\textwidth]{img/parameter_selection.pdf}
		\caption{Results of the PV selection by parameter removal for the RF model. The bars indicate the average MSE of all 100 RF with the PVs indicated. The line indicates the change in MSE between two bars. A bend in the line represents a strong decrease in the model accuracy which indicates that a PV crucial for the model was removed. For this example, the PVs GG, UP and YT were selected (compare table \ref{tab:VariableSelection}). Abbreviations: EA = electricity access, EC = energy consumption,  GG = greenhouse gas, UP = urban population, YP = young population and MR = MSW\textsubscript{recycled}.}
		\label{fig:VariableSelection}
	\end{figure*}
	\begin{table*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\rowcolors{17}{}{color2!10}
		\caption{Results of the PV selection. The numbers indicate the rank of the PV in the certain process (e.g. urban population was always the most important PV). A lower rank indicates a higher importance compared to high-ranked PVs. The PVs in bold were selected for the final model (compare Fig. \ref{fig:VariableSelection}).}
		\begin{tabularx}{\textwidth}{p{3cm}  >{\centering\arraybackslash}p{1.2cm} >{\centering\arraybackslash}p{1.2cm} c >{\centering\arraybackslash}p{1.2cm} >{\centering\arraybackslash}p{1.2cm} c >{\centering\arraybackslash}p{1.2cm} >{\centering\arraybackslash}p{1.2cm} c >{\centering\arraybackslash}p{1.2cm} >{\centering\arraybackslash}X}
			\toprule
			\multirow{2}{1cm}{\textbf{PV}}& \multicolumn{2}{@{}c@{}}{\textbf{Parameter removal}} & & \multicolumn{2}{@{}c@{}}{\textbf{Parameter shuffling}} & & \multicolumn{2}{@{}c@{}}{\textbf{Univariate}} & & Average & Median\\ \cline{2-3}\cline{5-6}\cline{8-9}
			& NNET & RF & & NNET & RF & & NNET & RF & & Rank & Rank \\ \toprule
			Electricity access & 4 & 4 & & 7 & 4 & & 5 & 7 & & 5.2 & 4.5\\
			Energy consumption & 5 & 5 & & 5 & 5 & & 1 & 3 & & 4.0 & 5\\
			GDP & 6 & 6 & & 3 & 7 & & 4 & 4 & & 5.0 & 5\\
			\textbf{Greenhouse gas} & \textbf{2} & \textbf{2} & & \textbf{2} & \textbf{2} & & \textbf{3} & \textbf{2} & & \textbf{2.2} & \textbf{2}\\
			\textbf{Urban population} & \textbf{1} & \textbf{1} & & \textbf{1} & \textbf{1} & & \textbf{2} & \textbf{1} & & \textbf{1.2} & \textbf{1}\\
			\textbf{Young population} & \textbf{3} & \textbf{3} & & \textbf{4} & \textbf{3} & & \textbf{6} & \textbf{5} & & \textbf{4.0} & \textbf{3.5}\\
			MSW recycled & 7 & 7 & & 6 & 6 & & 7 & 6 & & 6.5 & 6.5\\ \bottomrule
		\end{tabularx}
		\label{tab:VariableSelection}
	\end{table*}

	\subsubsection{Testing the model stability by varying the training data}
	\label{sec:MLMStability}
	After optimizing the model architectures and selecting the PVs, we were interested in the stability of the models with respect to the size of the training data and biases in the training data. Therefore, we separately trained 100 models using all the available data from each continent, cluster or GNI classification as biased training data. We selected the model with the lowest MSE to compare the effect of the selected bias in the training data. Additionally we compared the results of the training data size from models using randomly selected 90\%,80\%,...,10\% of the available data. These tests allow us to assess the robustness of the MLM to uncertainty and loss of the underlying data as well as the dependence of the MLM results on PV data substructures. The results of the stability tests are presented in section \ref{sec:ModelStability}.
	
	\subsubsection{Selection of the best model}
	\label{sec:SelectionOfBestModel}
	After optimizing the model architectures and selecting the PVs, we wanted to select the NNET and RF model that yields the best results. We therefore trained 10'000 NNET and RF models on randomly selected 90\% of the available data. We then used several techniques to assess the prediction performance, which will now be presented:\\
	\textbf{CV performance:} In a first step, each generated model is retained or excluded based on their performance in the CV. The CV was done while creating the model (section \ref{sec:MLMOptimization}). We only retain models whose MSE of the prediction using the unseen test-dataset is smaller than the MSE of at least one prediction using one of the nine training-datasets of the CV. This drops models that were not able to predict unseen data within the same MSE range as the training data used in the CV.\\
	\textbf{MSE:} In a second step, we rank the models by their MSE of the prediction of the test data. This selects models that were better in predicting the unseen data.\\
	\textbf{Errors for selected countries:} In a third step we ranked the models according to their prediction performance for selected countries. This is done to avoid large errors for countries with a big global impact. We ordered the models by their performance on the 20 countries that contribute the most to P\textsubscript{O} according to Jambeck et al. \cite{Jambeck2015}, and for which we had data available. This includes the eleven countries China, Indonesia, the Philippines, Sri Lanka, Egypt, Malaysia, Algeria, Turkey, Brazil, Morocco and the United States of America. The model prediction error of these countries was multiplied with the country's marine plastic input according to Jambeck et al. \cite{Jambeck2015}. This prefers models that show low prediction errors for the countries with the highest influence on P\textsubscript{O}.\\
	We then selected one NNET and RF model each. The selected model has one of the smallest MSE and its error for the prediction for the top 20 countries is also among the smallest. Additionally, its overall prediction error for all countries is as close to the smallest prediction error of all models. The selected NNET and RF models were finally used to predict the MSW\textsubscript{C} for each country and year. We present these models in section \ref{sec:PerformanceOfTheMLM}.\\
	
	The above described NNET and RF models are used to predict the MSW\textsubscript{C} for each of our 173 countries for the years 1990-2015. In order to create geospatial estimates in the WFM, the MSW\textsubscript{C} data has to be rasterized. This was done as described in section \ref{sec:Rasterization}. The procedure to estimate P\textsubscript{E} and P\textsubscript{O} from MSW\textsubscript{C} are described in the next section.
	
\subsection{Waste Flow Model}
\label{sec:WFM}
By running the MLM we estimated MSW\textsubscript{C} for each country in the time period 1990-2015. Next, we describe how we calculated how much plastic enters the ocean (P\textsubscript{O}) on a 0.1x0.1° raster using our WFM. First we present the calculation of the total amount of produced MSW (MSW\textsubscript{T}), consisting of the MSW inside the waste collection system (MSW\textsubscript{C}) and uncollected waste. Using this information, we then calculated the total amount of plastic waste (P\textsubscript{T}), again inside and outside the waste collection system. These steps are described in section \ref{sec:CalculatingMSWT}. Then, we show how we estimated the amount of plastic that is entering the environment (P\textsubscript{E}) in section \ref{sec:CalculatingPlasticEnteringEnvironment}. This includes both plastic lost from the waste collection system and littering as well as from uncollected plastic waste. Finally, we present how we estimated the amount of P\textsubscript{E} that enters the ocean by using a distance based transport probability (section \ref{sec:ProbabilityOfTransport}). This process is also visualized in Figure \ref{fig:WFM}.
			
	\begin{figure*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\centering
		\includegraphics[width=\textwidth]{img/Overwiev.pdf}
		\caption{Graphical representation of the waste flow model (WFM). We use machine learning models (MLM) to predict the amount of municipal solid waste collected at the households (MSW\textsubscript{C}). This process is described in section \ref{sec:ModellingMSWProduction}. Using the fraction of population with access to the waste collection system (F\textsubscript{collection}) we estimate the amount of uncollected MSW. The sum of MSW\textsubscript{C} and uncollected MSW represents the total amount of MSW (MSW\textsubscript{T}) that is produced. This is described in section \ref{sec:CalculatingMSWT}. We then use the fraction of plastic (F\textsubscript{plastic}) to estimate the amount of plastic inside and outside the waste collection system (P\textsubscript{T}). Then, we predict the losses of P\textsubscript{T} to the environment (P\textsubscript{E}). The losses are dependant on each countries waste collection system as described by the fractions of waste treated by recycling (F\textsubscript{recycled}), in landfills (F\textsubscript{landfills}) or in dumps (F\textsubscript{dumps}). A more detailed description can be found in section \ref{sec:CalculatingPlasticEnteringEnvironment}. Finally, we estimate the amount of P\textsubscript{E} that reaches the ocean (P\textsubscript{O}). This is done by applying distance based probabilities as described in section \ref{sec:ProbabilityOfTransport}.}
		\label{fig:WFM}
	\end{figure*}
	\begin{figure*}[htbp]
		\renewcommand{\arraystretch}{1.2}
		\centering
		\includegraphics[width=\textwidth]{img/waste_flow.pdf}
		\caption{Waste flow models for high-income (A), upper middle-income (B), lower middle-income (C) and low-income (D) countries. The fractions of collected and uncollected plastic waste are the average of MSW\textsubscript{collection} in 2015 and the country specific values were considered in the model. The fractions of leakage in B and C represent values of China and the Philippines respectively and are taken from the OceanConservancy and the McKinsey Center for Business and Environment \cite{McKinseyCenter2015}. Fractions of leakage for A and D are assumptions. In A, on average 2.2\% of all produced plastic waste leaks into the environment while in B, C and D averages of 6.2\%, 22.7\% and 38.0\% are leaked respectively.}
		\label{fig:WasteFlow}
\end{figure*}

	\subsubsection{Calculating the total amount of MSW and the total amount of plastic waste}
	\label{sec:CalculatingMSWT}
	MSW\textsubscript{C}, as calculated by the MLM using the PVs, only includes the MSW that was collected in each year. It does not include MSW that was not collected at the households (e.g. due to a missing collection system or due to missing regulation on MSW treatment \cite{Hoornweg2012}). To get an estimate of how much waste is uncollected, we used the fraction of population with access to the waste collection system (F\textsubscript{collection}). The total amount of MSW (MSW\textsubscript{T}) was thus obtained by adding the uncollected waste (outside waste treatment system) to MSW\textsubscript{C} (inside waste treatment system). By further applying each country's fraction of plastic in MSW (F\textsubscript{plastic}) we calculated the total amount of plastic waste produced (P\textsubscript{T}) on a 0.1x0.1° raster (both inside and outside the waste treatment system). In the next section we describe how this was used to estimate P\textsubscript{E}.
	
	\subsubsection{Calculating the amount of plastic entering the natural environment}
	\label{sec:CalculatingPlasticEnteringEnvironment}
	To estimate how much of the plastic inside and outside the waste collection system enters the environment, we needed to account for the MSW treatment situations of each country. As represented in Figure \ref{fig:WFM}, collected plastic waste can either be recycled, incinerated or it can be deposited in both landfills (government operated and controlled) or dumps (sites that are not under government control) \cite{Hoornweg2012}. Uncollected plastic on the other hand can be burned locally or put into dumps \cite{Hoornweg2012}. Losses from uncollected plastic to the environment are more likely to occur and a study in China, Indonesia, the Philippines, Thailand and Vietnam concluded that approximately 75\% of plastic in the environment comes from uncollected plastic waste \cite{McKinseyCenter2015} while the remaining 25\% leaked from the waste collection system. We therefore had to treat collected and uncollected waste differently to estimate the losses to the environment. We assume that the flux from P\textsubscript{T} to P\textsubscript{E} from both inside and outside the waste treatment system is dependant on the countries income \cite{Hoornweg2012}. Thus, we treated each income class differently, which is represented in Figure \ref{fig:WasteFlow} A-D. Following we present how the fluxes from P\textsubscript{T} to P\textsubscript{E} were estimated for plastic inside and outside the waste treatment system.\\
	Unfortunately it was not possible to find precise information on the loss of collected and uncollected plastic. However, the study by OceanConservancy and the McKinsey Center for Business and Environment \cite{McKinseyCenter2015} investigated plastic leakage fractions in the Philippines and China in detail. They found that 84\% of plastic that enters the environment in China comes from uncollected MSW (Fig. \ref{fig:WasteFlow} B) while 74\% of leaked plastic in the Philippines comes from waste inside the collection system (Fig. \ref{fig:WasteFlow} C). To adapt these values to the other countries, each country was grouped based on its GNI (according to the world bank definition \cite{Hoornweg2012}): high-income (Fig. \ref{fig:WasteFlow} A), upper middle-income (Fig. \ref{fig:WasteFlow} B), lower middle-income (Fig. \ref{fig:WasteFlow} C) and low-income (Fig. \ref{fig:WasteFlow} D). The loss fraction of collected and uncollected plastic of lower middle-income countries was then calculated by multiplying the reported fraction of loss of the Philippines (17\% for collected plastic, 31\% for uncollected plastic ,Fig. \ref{fig:WasteFlow} C)) with a scaling factor. This scaling factor was derived by dividing the fraction of dumps and landfills of the country by the fraction of dumps and landfills in the Philippines (compare table \ref{tab:DataCompleteness}). Thus, lower middle-income countries that dispose more of their waste in landfills and dumps than the Philippines are assumed to have a loss fraction that is bigger than the Philippines. Upper middle-income countries were compared to China in the same way. However, for high and low-income countries, no data could be obtained. Thus, the assumption was made that high-income countries loose 2\% of their collected waste (littering rate used in Jambeck et al. \cite{Jambeck2015}) and 10\% of their uncollected waste (Fig. \ref{fig:WasteFlow} A), while low-income countries loose 30\% and 45\% respectively (Fig. \ref{fig:WasteFlow} D).\\
	With the estimates of P\textsubscript{E} we then predicted P\textsubscript{O} using a distance based probability. This is described in the next section.
	
	\subsubsection{Probability of plastic waste transport}
	\label{sec:ProbabilityOfTransport}
	In order to estimate P\textsubscript{O} based on P\textsubscript{E}, we applied a distance based probability of transport. This includes riverine transport as well as a land-based transport to the shores. These two types of transport will be explained in this section.\\
	
	\textbf{Land based transport}\\
	A raster containing the minimum distance of each land cell to the closest ocean cell was calculated using the \textit{distance} function from the \textit{raster} package \cite{Hijmans2017}. The maximum travel distance was set to 100 km and the probability of mismanaged plastic waste to reach the ocean was set to 0\% outside this range. Cells inside the range obtained a probability that rises proportionally to the distance from 0\% at the maximum distance to 100\% at the shore ($p_{transport} = 1 - s/100 km$ where $s$ is the distance to the shore in km). This approach is similar to Jambeck et al. \cite{Jambeck2015}, who included the mismanaged waste produced by the population within 50 km from the shore. They then applied conversion rates of 15\%, 25\% and 40\% from mismanaged waste to marine litter to estimate the possible range of marine input. In our approach, we calculate the probability to enter the environment by losses to the environment from MSW\textsubscript{T} (from both the collected and uncollected waste), as opposed to Jambeck et al. approach to assess conversion rates. We argue that plastic in the environment that is closer to the shore has a higher probability to enter the ocean than waste further inland, while Jambeck et al. assume that all of the mismanaged waste reaches the ocean. We expect our results to be similar to Jambeck et al. because the integrated likelihood corresponds to the transport of 100\% of the mismanaged MSW within 50 km from the shore ($\int_{0}^{100}\frac{x}{100}dx = \int_{0}^{50}1 dx$).\\

	\textbf{Riverine transport}\\
	Similarly to the land-based transport, plastic waste in a distance of 20 km is assumed to have a chance to enter the river and this chance is proportional to the distance to the river. Furthermore, waste that enters the river upstream is considered to have a smaller chance to be transported all the way to the ocean than waste that enters further downstream. This probability was set to 100\% at the river mouth and sinks gradually to 0\% at a distance of 1'000 km upstream ($p_{transport} = 1 - s/1000 km$ where $s$ is the distance to the ocean in km). This approach is similar to Jambeck et al. The probability of mismanaged waste of a certain cell to reach the ocean by riverine transport was thus calculated by multiplying the probability to enter the river by the probability to be transported by the river. Thus, mismanaged plastic waste produced near a river in countries that are far from the shore will have a lower chance to end up in the ocean compared to mismanaged plastic waste produced near a river in a country closer to the ocean.\\
	With our river transport model, all types of plastic that enter the river will be transported. In reality, not all plastics will float. Due to their densities, only PE-HD, PE-MD (both 0.94-0.97-0.92 g/ml \cite{Moret-Ferguson2010}), PE-LD, PE-LLD (both 0.89-0.93 g/ml \cite{Moret-Ferguson2010}) and PP (0.85-0.92 g/ml \cite{Moret-Ferguson2010}) will float in fresh water (1 g/ml) or salt water (1.025 g/ml \cite{Andrady2011}). PVC (both 1.16-1.41 g/ml \cite{Moret-Ferguson2010}) and PET (both 1.38-1.41 g/ml \cite{Moret-Ferguson2010}) are too dense. However, also these plastic types could float with a respective shape (e.g. PET bottles). As plastic products are often mixed with fillers which changes their density \cite{Andrady2011}, predictions of the floating behaviour of plastic waste are difficult. However, non floating plastics can still get transported. Hurley \cite{Hurley2018} found that flooding cleared 70\% of all microplastics in river sediments. Thus, we assume that all plastics get transported and do not differentiate between the types in riverine transport. By decreasing the probability of transport with increasing distance from the sea we include losses on the surface (dams, beaching, entrapment \cite{Kooi2018}) and sediments.\\
	
\subsection[Estimates of P\textsubscript{O} for Future Projections \newline and Policy Chance Scenarios]{Estimates of P\textsubscript{O} for Future Projections and Policy Chance Scenarios}
\label{sec:FutureProjections}
With the methods described so far in this chapter we can estimate the worldwide P\textsubscript{O} on a 0.1x0.1° raster for the years 1990-2015. Besides the evolution of P\textsubscript{O} and its geospatial variances, we also want to use the MLM and WFM to predict marine plastic debris inputs for the future and for several scenarios. These include future projections for the years 2030 and 2050 as well as policy change scenarios affecting the losses from uncollected waste and the waste collection system. We created the following scenarios:
\begin{itemize}[noitemsep,topsep=0pt]
	\item \textbf{Future projections for 2030 and 2050}: We created future estimates to assess the change of predicted MSW\textsubscript{C} for 2030 and 2050 compared to 2015, assuming business as usual. We used the UN population estimates from 2015 to 2050 \cite{UnitedNations2018Pop} to calculate annual growth rates of the population for each country. We then applied these growth rates to the gridded population density data of 2015 \cite{SEDAC2017}. Thus, we assume that the population of each grid cell grows equally and neglect urbanization. The PVs were extrapolated to the years 2030 and 2050 using linear fits of the available data from 1990-2015 (as described in section \ref{sec:MissingData}). Thus, we create the future projections with the MLM while the WFM represents 2015 conditions.
	\item \textbf{100\% collection rate for MSW}: We estimate the effectiveness of policy measures that aim to collect all MSW (MSW\textsubscript{C} = MSW\textsubscript{T}). Therefore, we run the MLM with the 2015 data but assume that all of the MSW is collected (F\textsubscript{collection} = 100\%) and no uncollected MSW is generated. Thus, we only change parameters of the WFM.
	\item \textbf{No losses in the waste collection system or no losses of uncollected MSW}: We estimate the effectiveness of policy measures that aim to avoid all losses from the waste collection and treatment system. Again, we run the MLM with all the 2015 data and assume, that all of the MSW\textsubscript{C} remains in the collection system and no loss to the environment occurs. Uncollected MSW still leaks to the environment. To compare this with the contrary, we additionally create a (unrealistic) scenario where none of the uncollected MSW enters the environment while losses from the waste collection system occur. Here we assume that all of the informally disposed MSW is discarded in an environmentally safe way. Again, this changes only parameters of the WFM.
	\item \textbf{Decreased transport distances}: We estimate the effects of higher efforts in the management of MSW\textsubscript{C} as well as higher rates of collection of waste that leaked to the environment (community cleaning) or removal of plastic waste from rivers. Therefore, we reduce the maximum distance over which plastic waste can reach the ocean by 50\% to 50 km from land to the ocean, 10 km from land to rivers and 500 km along rivers to the ocean. This also changes only WFM parameters.
\end{itemize}

We estimated P\textsubscript{E} and P\textsubscript{O} using these future projections and policy measure scenarios. The results as well as a comparison to the year 2015 can be found in section \ref{sec:ScenarioResults}. Next we describe how we predicted the advection of P\textsubscript{O} in the ocean using a Lagrangian simulator.

%\subsection{Sensitivity Analysis of Output}
%\label{sec:SensitivityAnalysis}
%We assessed the sensitivity of the model output to changes in the parameter values to see how global development and change can influence the plastic input to the ocean.

\subsection{Simulation of Plastic Particle Distribution by \newline Ocean Currents}
\label{sec:SimulationofPlasticDistribution}
In this section we describe how we estimated the distribution of plastic particles entering the ocean on a 0.1x0.1° grid (P\textsubscript{O}). This was done by running PARCELS (Probably A Really Computationally Efficient Lagrangian Simulator) \cite{Lange2017} using the current data described in section \ref{sec:DataSources}. The PARCELS project is an open source set of python classes to create particle tracking simulations using hydrodynamic fields. We begin by describing how we modified the raster containing P\textsubscript{O} to allow usage in PARCELS. This is followed by the description of the simulation set-up of PARCELS. In the end of this section we present how we analysed the PARCELS output.
	
	\subsubsection{Modification of the particle input location}
	A modification of our P\textsubscript{O} data (0.1x0.1° raster) was needed as the current data is available on a 0.25x0.25° resolution (compare table \ref{tab:GeospatialData}). This is a necessary step because the fine 0.1x0.1° raster includes bays that are considered as land in the coarser 0.25x0.25° current data. Particles released at raster cells of the 0.1x0.1° raster might not be inside a cell of the 0.25x0.25° current data. PARCELS would thus not advect these particles as they are considered to be on land. As the current data is also known to be inaccurate close to the shore \cite{GlobCurrentConsortium2015} we decided to release particles from cells that are adjacent to shoreline cells (e.g. further into the ocean). We thus created a 1x1° raster that contains the P\textsubscript{O} of all 0.1x0.1° raster cells. In this new raster, the plastic input of each cell is moved to the closest adjacent cell further away from the shoreline. We obtained the coordinates of each cell that represents at least 100 t of plastic waste from P\textsubscript{O}. These coordinates were then used in simulation. Next we describe how the simulation was set-up.
	
	\subsubsection{Simulation set-up}
	%We run the PARCELS simulation from 1993 until 2015. Unfortunately, current data for the years 1990-1992 is not available at the desired tempo-spatial resolution. However, in order to create somewhat realistic initial conditions, each simulation included a pre-run phase where the years 1990-1992 were simulated using current data from 1993 but the included particle loads and locations were as calculated in the models. This ensures that the simulation does not start from an "empty" ocean.\\
	We run the PARCELS simulation from 1993 until 2015. We therefore use the current data (only from 1993-2017) as described in table \ref{tab:GeospatialData}) and the coordinates of cells where at least 100 t of plastic enter the ocean (previous section). We add particles every 10 days according to the mass of plastic entering the ocean (100 t = 1 particle, 200 t = 2 particles, ...). Each particle is then advected every 10 minutes for one day using 4\textsuperscript{th} order Runge-Kutta and the new position is determined. A Brownian motion kernel is included to simulate turbulent transport as well as micro scale advection not represented in the 0.25x0.25° hydrodynamic field. Therefore, each particles position is changed by a random value between 0 and 100 m/day in both longitudinal and latitudinal direction. This also ensures that particles with the same start point would have slightly different paths to represent both sub-grid scale motions as well as differences in particle size, shape and buoyancy. At the end of each month, the simulated particle density in the 0.25x0.25° ocean grid is stored. Furthermore, every particle is removed from the simulation after three years. This sink represents the settling of fragments as described by Koelmans et al. \cite{Koelmans2017}. Because the model starts with an uncontaminated ocean in 1993 and it takes three years to reach an equilibrium (when the initial particles get deleted), the PARCELS simulation is able to produce meaningful results for the time period 1996-2015.\\
	In the next section we describe how we interpreted the results from the PARCELS simulation using the set-up described above.
	
	\subsubsection{Analysis of particle densities}
	As mentioned above, a particle in the PARCELS simulation represent 100 tons of plastic. This particle size was chosen to reduce the runtime of the simulation. Because we included a Brownian motion kernel, we assume that an increased number of particles (e.g. 1 particle representing 10 t of plastic) would create smoother density distributions. Therefore, we smoothed the densities obtained from PARCELS. We used the \textit{kde2d} function from the \textit{MASS} package \cite{Ripley2018} to derive kernel density estimates (KDE). KDEs assume that the probability of each particles location can be described as normally distributed, which represents the uncertainty of the Brownian motion. The sum of all particle probabilities represents the cumulative density function. Thus, the product of the KDE and mass of plastic represented by each particle (100 t) equals the amount of plastic in the ocean. We obtained the standard deviation of each particles normally distributed location using the \textit{bandwidth.nrd} function from the \textit{MASS} package \cite{Ripley2018}. This function uses the \textit{Gaussian approximation} to obtain a optimal smoothness of the prediction compared to the data. The obtained KDE was then multiplied by the amount of plastic that each particle represents to get an estimate of the plastic density in the ocean.\\
	
In the following chapter we will present the results of the analyses that were described in this chapter.