A population-based simultaneous fugacity model design for polychlorinated biphenyls (PCBs) transport in an aquatic system

Graphical abstract PCBs Exchange between Organism and Environment (New Design).

Where M i represents the mass of PCBs accumulated in compartment i; V i m 3 À Á is the volume of compartment i; Z i (mol/Pa m 3 ) is the fugacity capacity of compartment i; f i ðPaÞ is the PCBs fugacity which represents the level of PCBs in compartment i. Thus, the dynamic change of PCBs in compartment i is estimated as: The formula is transformed through partial difference to become, As shown in formula (3), the change of PCBs fugacity in compartment i has three general contributions: the PCBs mass variation (dM i =dt), the change in fugacity capacity (dZ i =dt), and the change in compartment volume (dV i =dt). All compartments related to PCBs transport involve at least one of these three general components. To determine the fugacity variation in the certain compartment, we need to separately define the process and parameters in each media according to their physical, chemical, and biological features. Furthermore, we need to define the exchanging terms among different compartments. Moreover, we need to apply a method to estimate the existing biomass/population volume in each biotic compartment for population scale study.

The fugacity capacity Air
The main components related to PCBs transport in the air are a gaseous phase and aerosol. The fugacity capacity in each component could be expressed as [2]: Where R is the ideal gas constant (8.314 J/mol K); K OA is the octanol-air partition coefficient; T A is the air temperature (K). Thus, the fugacity capacity in air compartment should be: Where t 1 is the volume fraction of aerosol.

Water
The water compartment contains water column and suspended sediment. In some model designs, water compartment also includes aquatic organisms. However, since organisms are isolated and calculated separately, we isolate organisms from the water compartment. The fugacity capacity could be expressed as: Where H is the Henry's Law constant (Pa m 3 /mol); r 5 is the suspended sediment density (kg=m 3 ); d 5 is the mass fraction of the organic carbon; K OC is the organic carbon partition coefficient (L/kg), which is approximately 0.41 times of the Kow [3]: Thus, the capacity of water compartment should be: Where t 2 is the volume fraction of suspended sediment.
According to assumption 1, the water phase and the suspended sediment particle should have identical fugacities during PCB transport. In 2011, LimnoTech provided a study report regarding PCB loading patterns in Lake Ontario [4]. According to the study, the PCB input of Lake Ontario in 2005 came from air transmission (20%) and water flows (80%). Moreover, a detailed analysis on aquatic PCBs input indicates a 70%/30% allocation between dissolved PCBs (water column) and particle PCBs (suspended sediment Thus t 2 ¼ 3000

Sediment
The vertical homogeneity conversion of the sediment compartment is difficult since the PCBs content with sediment depth depends on the PCB contamination level during the deposition period. As a temporary compromise, the current sediment compartment only includes the very top layer of bio-active sediment (%0.1 m) which contains 10% dry residual mixed with around 90% of the saturated water (in volume fraction). The sediment is considered as flooded sediment, where little air existed in the compartment [19]. As a result, the sediment compartment is a mixture of water and sediment solid with organic particle attached to the organic matters. The solid sediment particle phase has a fugacity capacity as: Thus, the fugacity capacity of the sediment compartment is: Where t 3 is the volume fraction of solid sediment.
For accurate estimation, the fraction of organic carbon in flooded sediment could be calculated through water content and dry bulk density [5]: Where OC is the organic carbon concentration (mg/dw g). The inorganic sediment particle density is conventionally taken 2.65 g/cm 3 ; the density of organic matters can be corrected assuming a density of 1.25 g/cm 3 . Thus, the sediment solid density can be expressed as: Thus, the water content is: Finally, the fraction of OC is: 1:776 À 0:363lnð1000ÃOM=1:7Þ 2:65 À 1:4% OM This equation means we can use water content to estimate the fraction of organic carbon in the sediment.

Organism
According to Mackay, the fugacity capacity of biota is defined as [2]: Where L is the lipid fraction in the organism.

PCB mass variation
The PCB mass variation, or dM i =dt, is defined as the PCB mass enters or exits the system with general transport processes. The general form for the changes of fugacity in compartment i can be expressed as, In this formula, i represents the different media; j represents other media that interact with media i; M i mol ð Þ represents the current PCB mass in medium i; tðdayÞ represents the PCB transport and allocation time; E i ðmol=dayÞ represents the direct pollution exchange rate to medium i; f j ðPaÞ represents the PCB fugacity in medium j; f i Pa ð Þ represents the fugacity of medium i; D ji (mol/Pa day) represents the PCB transport processes from medium j to medium i ði 6 ¼ jÞ; D T i (mol/Pa day) represents the total PCB elimination/exit from medium i.

Air
Pollutant transport processes related to air compartment include three processes: the inter-media exchange, the self-elimination, the systematic exchange [6]. During the inter-media exchange, the entrée is mainly through water volatilization (air-water diffusion, D V ), while the exit pathways include absorption (water-air diffusion, D V ), wet dissolution (D RWW ), dry deposition (D QDW ), wet particle deposition (D QWW ). Since no biota is considered in the air, no direct exchange between the air compartment and any organisms. The self-elimination, or reaction (R A ) within the compartment generally eliminate contaminate through photodegradation and is related to the compartmentalbased lifetime. Finally, the systematic exchange is mainly through the advections (D AI =D AO ). As a result, the fugacity variation in the air compartment could be written as: Where Diffusion : Wet Dissolution : Dry Deposition : Wet Particle Deposition : Reaction : Advection Input : Advection Output : Thus, Biota The organisms are discussed previously for better understanding their interactions with the environment phases. In this study, the definition of the inter-exchange process among different biota groups occurs only within the food web, while the processes with the environmental groups are identified as a systematic exchange. In the inter-exchange process, PCBs are absorbed by organisms through food ingestion (D FI ) and are released through predation (D Pred ).
When studying the PCB transport between environment and organism in water compartment, organisms are divided into pelagic and benthic species, because habitat location will lead to different calculation method PCB exchange rate. Gill uptake is one of the primary routes to transfer PCBs into organisms (D GGW =D GGS ) [18]. The pathways to transport PCBs to the environment includes gill release (D GLW =D GLS ), natural mortality (D MD ), and egestion (D E ). Egestion is combined by the undigested food (1 À E D ) and PCB exchange between the gut and the feces (D Ethe X ). Undigested food is usually estimated as a proportion of the total food ingestion, while the gut/feces exchange rate is estimated through trophic magnification factor (TMF) and trophic levels [16]. The PCB self-elimination in biota group is mainly through metabolism (R B ).
Considering the existence of the decomposing process, we assume that the PCB inside dead organisms caused by natural mortality will be initially decomposed and released to the environment before regaining through the food web. For the pelagic species, PCBs from decomposed organisms return to both water and sediment; for benthic groups, all released PCBs go to the sediment compartment. Thus, the changes of fugacity in biota could be expressed as: Pelagic species Benthic species Where Gill Uptake : Food Ingestion : Gill Release : Metabiolism : PCB exchange between the gut and the feces The PCB exchange rate between the gut and the feces can be calculated through TMF and trophic levels. TMF, or trophic magnification factor, could be used to evaluate the proportion of PCB escape from the system through feces [16]. The species trophic level can be calculated by the following formula [7]: Where TL j represents the fractional trophic level of prey j, and p ij represents the fraction of j in the diet of i. The PCB released through the fence is then decided by the true TMF differences between food and diet: Predation.
Natural mortality (mortality without predation) The natural mortality rate is estimated by Then et al. in 2015, who used over 200 fish species to evaluate the current existing empirical models for natural mortality rate estimation [8]. We selected one of the best models as a basis for estimating the natural mortality loss.
Where, t max is the maximum surviving time for species i (years);

Growth dilution
According to the new representation in formula (1) through (3), the growth dilution does not belong to the first category, dM i dt , since there is no actual entrée or exit of any PCB during the process. It is merely a volume change. As a result, it should be moved to the third part.
In sum, the extended expressions for formula (29) and (30) are: Pelagic species Benthic species Water The pollution transport through water is more complex than the air section because of the existence of organisms. To achieve fidelity, PCB exchange processes between the environment and organisms have been built. Similarly, the pollutant exchange in the water section is divided into three parts. The entrée processes in intermedia exchange include The self-elimination, or reaction (R W ) within the compartment generally eliminate contaminate under a first-order decay rate, which is relative to its compartment-based lifetime. Finally, the systematic exchange is mainly through the advections in/out (D WI =D WO ) of the system. As a result, the fugacity variation in the water compartment could be written as: Advection in=out : Water À sediment Dif f usion : Water À Sediment Deposition : Water À Sediment Resuspension : Reaction : In sum the extended expression for water compartment could be written as: Sediment Similarly, as water compartment, the sediment compartment includes biotic activities. Thus it also contains similar processes. The entrée processes in intermedia exchange include The self-elimination, or reaction (R W ) within the compartment eliminates PCB through biodegradation (aerobic remediation only), which is relative to its compartment-based lifetime. The sediment compartment does not have a direct PCB input route, but sediment compartment can push PCBs out of the system by deposition. As a result, the fugacity variation in the sediment compartment could be written as: Deposition : Reaction : In sum the extended expression for sediment compartment could be written as: Fugacity capacity variation

Air
The fugacity capacity of air compartment could be expressed as: According to Li et al., the Octanol/Air partition coefficient is temperature sensitive with an estimation of: Thus, the fugacity capacity variation in air compartment is expressed as: For Water The fugacity capacity in water could be expressed as: Where H is Henry's law constant. According to research by Schwarzenbach in 2002, the Henry's Law constant could be affected by the temperature with the following formula, also as known as van't Hoff correction [9,10]: Where H ref is the referenced Henry's Law constant at T ref ; U AW is the is the difference in internal energies of PCB in the phase change from air to water (kJ=mol). Similarly, the K OW also adept in van't Hoff correction: Where DU OW is the internal energies requirement for PCB going from octanol to water (kJ=mol).
However, the K OW is much less sensitive to the temperature variation. In this study, we can assume a constant K OW to simplify the calculation. Thus, the fugacity rate of change in water is: And Finally,

Sediment
The fugacity capacity in sediment could be expressed as: Thus, Where r 4 is the sediment density (kg/L). Thus the fugacity capacity change in sediment is: Finally, Biota According to Mackay, the fugacity capacity of biota is defined as: Where L is the lipid fraction in biota, then Then Compartment volume variation

Environmental compartment
In this study, we assume no volume change occurs in the environmental compartment.
Organism volume estimation through natural mortality, growth rate, and predation The organism population size is another essential factor for the new model. Pre-existing methods for biomass size estimation involve field investigation and measurement. To improve the estimation of biomass volume, we develop an energy-mass method, which could estimate the primary producers' biomass in most trophic structures. In this method, we use energy flow to estimate the biomass of the primary producers and mass balance to estimate all following species.
The estimation begins with primary producers. To begin with, two assumptions must be considered. First, it is assumed that the plant cover, such as grassland, forest, wetland, and so on, has long been existing and stabilized before the estimation. Second, solar energy is identified as the primary energy source in the ecosystem. Solar energy input is the fundamental source to determine the scale of the primary producer. However, it is not the only dominant factor to control the plant biomass. The amount of the plant biomass is also determined by the nutritional condition, water supply, and for aquatic ecosystems, light penetration. In general, most of these parameters cannot be measured directly, but they can be estimated via other parameters, such as photosynthetic efficiency and vegetation coverage [11]. Based on the pre-assumptions, the following formula is used for plant scale estimation: Where E Solar (J/s m 2 ) represents the total solar energy input to the unit surface; ' PE i % ð Þ represents the photosynthetic efficiency of plant i; s C i % ð Þ represents vegetation coverage rate of each type of plant in the study area; ' T i % ð Þ is the energy transport factor, the efficient proportion of energy stored in the system; # C i ðg=JÞ represents the carbon production factor which is the energy transferred to carbon in the system; ' C i ð%Þ represents the carbon fraction, that is, the weight percentage of carbon in the target organism i; t i days ð Þ represents the average lifetime of the species; A(m 2 ) is the surface area. In formula (69), the total solar radiation is acquired from the Solargis [12]. Notice that in the current study, the Photosynthetic Efficiency and Carbon Production Factor are usually measured together. The combined parameter of Photosynthetic Efficiency and Carbon Production Factor is based on the Green Solar Collector; converting sunlight into algal biomass [13]. The estimation of the biomass growth efficiency not only involves the photosynthesis efficiency but also takes into account daily consumption for organism growth and self-maintenance. The vegetation coverage rate can also be found in books [14] and the USGS GAP Land Cover Data Set. The calculation results are expressed as volume or mass since the density of most aquatic organisms is near water density. According to the observation data in Lake Ontario [15], the biomass density of phytoplankton is around 0.01-1 g/m 3 . The formula (69) calculation results, depending on the coverage rate and seasonal features, are around 0.03-1.5 g/m 3 . The next step is to calculate the higher trophic level biomass. Since all consumers gain their energy through food ingestion, it is continent to use food mass flows to monitor the biomass in each species. In current model design, we assume constant population sizes among all the biotic compartments in the ecosystem. The population size could increase through growth/reproduction and lose its size through natural mortality and predation. We do not consider disasters or incidents which could dramatically alter the population scale. Thus, the food mass flow could be written as: Where G represents growth rate (day À1 ) [18], P represent predation rate (day À1 ), and M is the mortality rate (day À1 ). The expressions for growth rate, predation rate, and natural mortality rate are: If we knew the food web details and the size of primary producer, it is easy to get all the rest populations as long as they are connected to the food web. In contrast, if the population of the species is known, we can also use this formula to calculate the correct food composition of each species.

Data validation and case study
We apply the new model to a series of case studies on polychlorinated biphenyl (PCBs) transport in Lake Ontario. The results are published in a recent article in Environmental Pollution: Modeling the impact of biota on polychlorinated biphenyls (PCBs) fate and transport in Lake Ontario using a population-based multi-compartment fugacity approach [20]. (https://doi.org/10.1016/j.envpol.2018.05.068).