An Underwater Light Attenuation Scheme for Marine Ecosystem Models 5 a

Simulation of underwater light is essential for modeling marine ecosystems. A new model of underwater light attenuation is presented and compared with previous models. In situ data collected in Monterey Bay, CA. during September 2006 are used for validation. It is demonstrated that while the new light model is computationally simple and efficient it maintains accuracy and flexibility. When this light model is incorporated into an ecosystem model, the correlation between modeled and observed coastal chlorophyll is improved over an eight-year time period. While the simulation of a deep chlorophyll maximum demonstrates the effect of the new model at depth. © 2008 Optical Society of America OCIS codes: (010.4450) Oceanic optics; (000.4430) Numerical approximation and analysis; (280.4991 Passive remote sensing References and links 1. B. Penta and J. J. Walsh, "A one-dimensional ecological model of summer oxygen distributions within the Chukchi Sea," Cont. Shelf Res. 15, 337-356 (1995). 2. O. Ziclinski, O. Llinas, A. Oschlics, and R. Rcutcr, "Underwater light field and its effect on a oncdimcnsional ecosystem model at station RSTOC, north of the Canary Islands," Deep Sea Res. 49, 3529-3542 (2002). 3. B. Penta, "Phytoplankton competition on the West Florida shelf: a simulation analysis with red tide implications," Ph.D. Dissertation, (Univ. So. Fla. 2000). 4. J. J. Walsh, B. Penta, D. A. Dictcrlc, and W. P. Bissctt, " Predictive ecological modeling of harmful algal blooms," Hum. Ecol. Risk Asses. 7, 1369-1383 (2001). 5. C. D. Moblcy, "Hydrolight 3.0 users' guide," (SRI Int., Mcnlo Park, Calif. 1995). 6. Z. Lee, K. Du, R. Arnonc, S. Licw, and B. Penta, "Penetration of solar radiation in the upper ocean: a numerical model for oceanic and coastal waters," J. Geophys. Res. 110, C09019 doi: 10.1029/2004JC002780 (2005). 7. B. Penta et al.. arc preparing a paper to be called "Ecosystem dynamics in the California Current System: analysis from a coupled physical-ecological ocean model." 8. I. Shulman, J. Kindle, P. Martin, S. dcRada, J. Doyle, B. Penta, S. Anderson, F. Chavez, J. Paduan, and S. Ramp," Modeling of upwclling/rclaxation events with the Navy Coastal Ocean Model," J. Geophys. Res. 112, C06023, doi: 10.1029/2006JC003946 (2007). 9. I. Shulman. J. Kindle, S. dcRada, S. Anderson, B. Penta, and P. Martin, "Development of a hierarchy of different resolution models for study U.S. West Coast California Current Ecosystem," in Estuarine and Coastal Modeling, M.L. Spaulding, cd. (Proceedings of the 8 International Conference on Estuarine and Coastal Modeling, 2004). 10. M. J. R. Fasham, H. W. Ducklow, and S. M. McKclvic, "A nitrogen-based model of phytoplankton dynamics in the oceanic mixed layer," J. Mar. Res. 48, 591-639 (1991). 11. F. Chai, R. C. Dugdalc, T.-H. Peng, F. P. Wilkcrson, and R. T. Barber, "One-dimensional ecosystem model of the equatorial Pacific upwclling system. Part I: model development and silicon and nitrogen cycle," Deep Sea Res. II 49, 2713-2745 (2002). 12. R. M. Hodur, J. Pullcn, J. Cummings, X. Hong, J. D. Doyle, P. J. Martin, and M. A. Rcnnick, "The Coupled Ocean/Atmosphere Mcsoscalc Prediction System (COAMPS)," Oceanography 15, 88-98 (2002). 13. P. A. Ncwburgcr, J. S. Allen, and Y. H. Spitz, "Analysis and comparison of three ecosystem models," J. Geophys. Res. 108, doi:IO.I029/200IJCOOI 182 (2003). #99635 $15.00 USD Received 31 Jul 2008; revised 10 Scp 2008; accepted 12 Scp 2008; published 2 Oct 2008 (C) 2008 OSA 13 October 2008 / Vol. 16, No. 21 / OPTICS EXPRESS 16581


Introduction
The underwater distribution and variability of photosynthetically active radiation (PAR) are important components in numerical models of marine ecosystems [1].Simulated primaryproduction and the resulting biomass distribution are dependent upon the numerical scheme used to describe the attenuation of PAR with depth [2], Variations of the euphotic zone depth can have large effects on the gross integrated primary productivity.Light limitation can influence the competition between phytoplankton groups [3,4] and the resulting species (or functional group) composition, as well as the vertical distribution of phytoplankton [2].
Spectral radiation models can reproduce the underwater light field with high accuracy [5], However, these models are computationally intensive and therefore, are inefficient for incorporation into large-scale, high-resolution, three-dimensional coupled models with short time-steps.An accurate, yet computationally simple model for light is needed.
Lee et al. [6] describe an innovative method for simulating the penetration of solar radiation in marine waters from remote sensing (ocean color) data.We have adapted this method for inclusion into a coupled biological-physical modeling system and validated its use with in situ data collected in Monterey Bay, CA during September 2006.

Models
The NCOM-CCS model (Navy Coastal Ocean Model of the California Current System) is a coupled bio-physical model for the West coast of the United States [7][8][9].The circulation model is based on Navy Coastal Ocean Model [8,9].The model ecosystem [7,10,11] consists of three state-variables for nutrients and two state-variables each for phytoplankton, zooplankton, and detritus.Phytoplankton photosynthesis is driven by PAR, which is derived from the high-resolution surface shortwave radiation fluxes from the Coupled Ocean-Atmosphere Mesoscale Prediction System (COAMPS•) [12].PAR at depth zero is computed as the fraction (0.48) of total solar irradiance that penetrates below the air-sea interface.Our preliminary treatment of PAR attenuation applied the Lambert-Beer law using attenuation coefficients for seawater [including chromophoric dissolved organic matter (CDOM), detritus, and total suspended solids (TSS)], k w , and chlorophyll k^: In this expression, all of the wavelengths that make up PAR (400 nm -700 nm) are attenuated equally with depth.Using Eq. ( 1), four different sets of attenuation coefficients (SI-4) were tested in the present study (Table 1).The values of SI are taken from the equatorial Pacific upwelling model of Chai et al. [11].S2 values were determined by Newberger et al. [13] for an upwelling region along the Oregon coast.Olivieri and Chavez [14] computed S3 from chlorophyll and attenuation measurements in Monterey Bay, CA.And Fujii et al. [15] used a multispectral optical model coupled to the ecosystem model of Chai et al. [11] to tune the attenuation coefficients of S4.The scheme of Lee et al. [6] was developed to utilize the inherent optical properties (IOP) absorption (a) and backscattering (b b ) of marine waters (at 490 nm) with the idea that these inherent optical properties can be retrieved via remote sensing.For total absorption at 490 nm in a range of 0.015 to 1.5 m" 1 (equivalent to chlorophyll concentrations ranging from 0.03 to 30 mg m 3 ), Lee et al. [6] (using HydroLight [5] simulations) developed a simple model for the vertical transmittance of solar radiation in the visible (i.e., PAR: 350-700 nm) domain, with a, bb, solar angle (8 a ), and depth (z) as variables: A key feature of this approach is that the attenuation of PAR is no longer treated as a vertical constant, but represents the change of light quality with increasing depth (i.e., longer wavelengths are attenuated rapidly in the surface water while the shorter wavelengths penetrate deeper).Lee et al. [6] assumed the IOP products derived via remote sensing were homogeneously distributed in the vertical and recommended caution when applying the satellite-based products to vertically stratified regions or depths below the well-mixed surface layer.Hereafter, this approach is called the "Rl" scheme.

Adaptation of the light scheme for non-homogeneous vertical distributions
We have extended the scheme of Eq. ( 2) to incorporate non-homogeneous vertical distributions of IOPs.Now, rather than extrapolating above surface observations to a homogeneous vertical distribution, the vertical profiles derived from observations or model predictions are used.This implementation allows the modeling of PAR attenuation in stratified regions and at depths below the ~ Me depth of satellite detection.
The absorption and backscattering coefficients in the new scheme (called Ml) are calculated from seawater (a w (490) = 0.015) [16] and the phytoplankton chlorophyll (observed or model-derived).Phytoplankton absorption is computed using the chlorophyll specific absorption coefficient at 490nm (a*(490) = 0.0375) [17]: where chl(z) are observed or model-predicted chlorophyll profiles.
For this first analysis we used an a* value based on field measured absorption coefficient and chlorophyll concentration and assumed that remains constant in the course of the modeling.In the future, an a* value that varies with phytoplankton functional group and/or with light and/or nutrient history could also be included when the necessary information is available [18].Scattering (b) is computed from this phytoplankton chlorophyll as well [19], then converted to backscattering by assuming b b as 1.5% of b [20].As backscattering coefficient makes minor contribution to the attenuation of PAR for waters in this study, errors in bb to b ratio have negligible impacts to the modeled PAR profiles.The solar angle is computed from latitude, longitude, date, and time of day.The spectral effects clouds have on PAR are minimal [21]; therefore, the model parameters developed in Lee et al. [6] are applicable to overcast sky conditions.Under these conditions, the sun angle is set to 45 degrees [22].

Data
Our in situ validation data set was compiled from a series of 98 profiles measured during 5-15 September 2006 in Monterey Bay, California (Fig. 1) aboard the R/V John H. Martin.Hyperspectral measurements were collected using a Satlantic HyperPro II profiling radiometer package equipped with a two-channel backscatter and fluorescence sensor (WETLabs ECO-BB2F; b b :470 and 700 nm, chlorophyll fluorescence).Temperature, conductivity (converted to salinity) and depth were recorded simultaneously by the profiler.A spectrally matched planar irradiance meter was used as a deck reference for all casts.A minimum of three up-and down-casts were collected at each station, with poor casts (ship shadow, high tilt angles, or excessively noisy data) discarded.
Data were processed using Prosoft 7.7.9(Satlantic, Inc.) with calibration and instrument files generated by Satlantic immediately prior to the field campaign in 2006.Data were processed (level 3) using shutter (dark) correction, depth binned to 1 m with 1 nm spectral resolution.Subsequent processing (level 4) utilized the following values.Integration Points: 5; Reflection Albedo: 0.043; Reflectance Index: 0.021; Refractive Index: 1.345, to produce the diffuse attenuation coefficient, downwelling irradiance, upwelling radiance, fluorescence, backscatter, water temperature, salinity, and density values used in this study.Processed Level 3 and 4 data were extracted to MATLAB (The Mathworks, Inc.) data format using the MAT Data Extractor in ProSoft.Remote sensing data, derived from the Moderate Resolution Imaging Spectroradiometer (MODIS) aboard the Earth Observing System (EOS) spacecraft Aqua, was processed using the Naval Research Laboratory's automated processing system (APS v5.4) [25].Chlorophyll derived from the OC3m algorithm [26] compare well with chlorophyll derived from the in situ fluorometer.The cloud-free mean field, for the time frame of the in situ observations, of (oc3m) chlorophyll is shown in Fig. 1(a).The other satellite derived products used in this study, computed at 488mn with the APS's quasi-analytical algorithm (QAA) [27]  During the time period of the observations, winds were predominantly from the North-Northwest creating upwelling favorable conditions [28].However, during 7-9 Sept., the winds switched to the South-Southwest.Also during this time period, the northeastern part of Monterey Bay experienced an unusually large phytoplankton bloom.The bloom consisted of Akashiwo sanguinea (Gymnodinium senguineum, G. splendens), a large (40-75 um), mixotrophic dinoflagellate species, known to form dense blooms in the Eastern Pacific [29].This species occurs in coastal waters offshore of San Francisco Bay during autumn or during periods of diminished upwelling when stratified conditions favor vertical migrators [29].The blooms can then be transported to San Francisco Bay [29].There were two relaxation events just prior to the time frame of these observations, 25-27 Aug. and 31 Aug. -2 Sept. The highest chlorophyll concentrations (>60 mg m 3 ) were measured (in situ) on 7 Sept. Upwelling favorable winds returned on 10 Sept. and the bloom began to dissipate between 12-15 Sept.These protozoans generally cannot migrate across strong density gradients.The vertical position of the dinoflagellate bloom [determined by their particle volume fraction (uL/L) to particle diameter (um) ratio] appears bounded by density gradients.

Methodology
To validate Ml, the vertical distributions of chlorophyll were estimated from in situ fluorometer observations described in the Data section.Using Eq. (3) these profiles were then converted to profiles of absorption (a) and backscattering (bb) at 490 nm.The upper most measured chlorophyll concentration was simply treated as the surface concentration.The sun angle computed from the time and position information was used except during overcast conditions (when 9 a is set to 45°).The light profiles obtained from Ml and the standard method (SI-4) were then compared to PAR measured concurrent with the fluorometer data.Next, using a(488) and b b (488) from the QAA processed data (or the temporal mean field [Figs.1(b)-1(c)] when the in situ profile was obscured by clouds), we ran the IOP-based scheme (Rl) for the time and location of each of the in situ profiles.
Lastly, we ran two Rl/Ml hybrid scenarios (M2-3).It is known that the absorption by CDOM is an important component in determining the attenuation of light [15,30].In the first scenario, we simply added the absorption due to CDOM and detritus, a<j g (488), derived from the QAA [Fig.1(d)] and applied it as a constant with depth as in the satellite model.In the second case, based on profiles of bb (Fig. 2) and b b :b ratios from the in situ data set, we forced the satellite derived CDOM (and detritus) absorption coefficients to have profiles matching the fluorometer profiles shape (i.e., we assumed case I like dependency).

Results
Each of the 98 in situ fluorometer profiles was used as input to estimate downwelling PAR with each of the eight light schemes (SI-4, Rl, and Ml-3).The outputs from the eight schemes were then compared to the in situ PAR profiles (as % of surface PAR).These comparisons are shown for two of the profiles, one with high chlorophyll and low PAR penetration [Fig.3(a)] and one with a more typical chlorophyll concentration profile [Fig.3(b)].Mean profiles were created from all the data and each of the models.Because the number of deep profiles is very small, the mean profiles are computed to 35 m depth only.Next, the error between each model mean % PAR profile and the mean measured profile are normalized by the data and the results plotted in Fig. 4(a).The mean fluorometer profile and number of stations with measurements taken at each depth are shown in Fig. 4(b); it is important to note the large drop-off of sample size around 30 m.The root mean square error (RMSE) and correlation coefficients of the various attenuation models (comparing the means) are given in Table 2.The shapes of the mean attenuation profiles are all generally well correlated with the mean data profile since they are all exponential decay functions.The M1&3 and S2-4 all match the shape of the mean profile with correlation coefficients > 0.99.Small offsets of the model profiles from the data may be due to the simple method used to extrapolate the fluorometer measurements to the surface.It is apparent from Fig. 4 that the standard (SI-4) model over-estimates PAR in the upper water column with each of the four sets of parameters.The longer (red) wavelengths of the PAR spectrum should be removed rapidly in this region.An important feature of the profiles from the standard model is switch-over at depth (the exact depth depends on the set of coefficients used) from the over-estimation of light to its under-estimation; now the remaining blue-green bands should not be absorbed so strongly.S4, with the coefficients optimized via the bio-optical model [15] has the smallest RMSE (0.63) of the standard formulations.The Rl model computed too-little light for all but the near surface because of its assumption of a homogeneous profile (i.e., continuous high chlorophyll for the whole water column).In the upper-most water column, in the region where the high chlorophyll values are approximately homogeneous, the Rl model simulates PAR well (Fig. 4).The RMSE of this model, 0.69, is comparable in magnitude with S4 but with too little rather than too much light.The MI model RMSE is 0.32 and it is apparent from Fig. 4 that the error is approximately constant with depth with no gross over or under-estimates.The inclusion of the remotely sensed CDOM in M2 as a homogeneous concentration with depth causes the modeled %PAR profile to behave like the Rl model profile with a smaller underestimation and a RMSE of 0.47.If, as with these data, the use of a Case 1 water like dependency is warranted, M3 provides the best results of all the models evaluated.M3 produces results with a small and non-depth dependent RMSE of 0.16.
In order to investigate results from the various light models at depths greater than those attained by our in situ sampling, we used a simulated 250 m water column.As evidenced by the straight lines when plotted on a semi-log plot (Fig. 5), the standard model contains no information about the spectral nature of the attenuation of light by water (the chlorophyll concentration was set to zero for the entire water column for these model runs).The IOP based models (represented by Ml) display a changing slope with depth, indicating the "pseudo-spectral" nature of the scheme.Remember, the values of the coefficient a w used in the standard models contain varying contributions from CDOM, detritus, and TSS depending on the coefficient used (see Models section).The dashed line marks the one-percent light level, used to define the base of the euphotic zone.While this figure may be an exaggeration of the effect, the differences in depth of the euphotic zone between models is apparent.

Effects of light schemes on the ecosystem model
The consequences of simulated euphotic zone depth variation between the light schemes (S2 and Ml) are evident with the incorporation of the schemes into the coupled ecosystem model.Phytoplankton utilize photons incident from any direction.Therefore, for completeness, downwelling PAR needs to be converted to scalar PAR for use in photosynthesis computations, which is achieved with division by the average cosine (/7).Average cosine, in general, changes within -0.5 to 1.0 vertically [23,24].In this study, for simplification, the constant value 0.7 is used for the average cosine.The small error in treatment of this parameter has no adverse impact on the overall conclusions because the same parameter is used in all simulations.
The NCOM-CCS simulation with the Ml light scheme produces a deep chlorophyll maximum (DCM) at the nutricline, which is not replicated in simulations using the standard  6).The NCOM-CCS simulated DCM is similar that observed along 140° W between 25° N and 35° N by Hodges and Rudnick [31] Jan. -Feb.1997.The ability to model this ubiquitous oceanographic feature [32] is an important consequence of our new scheme.Near-coast (100 km), monthly mean surface chlorophyll concentrations from NCOM-CCS (with both S2 and Ml light schemes) from 30° N to 46° N were compared to Sea-WiFS derived values over eight-year simulations (1999 -2006).The correlation between the model and satellite derived values is enhanced by the use of the new light scheme (Fig. 7).

Conclusions
The new scheme of underwater PAR propagation presented here, and in Lee el al. [6], has several advantages for use in coupled ecosystem models.This new scheme (Ml) is both accurate and expeditious.It is also flexible; attenuating and scattering components can be added either within the ecosystem model, or from remote sensing (M2), or other data sources (M3) to improve the accuracy of the modeled PAR field.The simulation of three-dimensional chlorophyll concentrations within the NCOM-CCS coupled model is improved (Figs. 6-7) using this new scheme.
Fig. I. Mean (of cloud-frcc pixels) OC'3m chlorophyll (a), total absorption (b), backscattering (c), and absorption due to CDOM and detritus (d) in Monterey Bay, CA. for the period September 5-15, 2006.Symbols indicate the locations and dates of in situ profiles.
are: total absorption [Fig.1(b)], backscattering [Fig.1(c)], and the absorption due to CDOM and detritus [Fig.1(d)].The dates and locations of the in situ hydrocasts are overlain on the fields in Fig 1.

Fig. 3 .
Fig. 3. PAR profiles (%) from models and in situ data with in situ fluorometer profile for two example casts.Low (a) PAR penetration/high chlorophyll and deeper (b) PAR penetration/lower chlorophyll.

Fig. 4 .
Fig. 4. Normalized % error of the model means compared with the data mean profile (a).Mean in situ fluoromctcr profile and total number of samples from each depth (b).

Fig. 5 .
Fig. 5. Comparison % PAR profiles for SI-4 and Ml run for an idealized water column with no chlorophyll.Dashed light delineates the 1% light level (i.e., depth of the euphotic zone).

Table 1 .
Coefficient sets used for the standard PAR attenuation model.

Table 2 .
Results of model/data comparisons.