Final_Project_Proposal.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Final Project Proposal\n",
    "\n",
    "Keith Arora-Williams<br>\n",
    "March 21, 2016<br>\n",
    "Landscape Hydrology"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Background"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Microbial Metabolism and Meteorological Process\n",
    "\n",
    "I am interested in how meteorological, seasonal, and physical processes affect microbial metabolism within a lake. A biogeochemical model developed by Hunter at al. <a name=\"ref-1\"/>[(1998)](#cite-Hunter199853) and modified by Prehiem et al. <a name=\"ref-2\"/>[(2016)](#cite-Preheim2016) was implemented to predict the metabolic redox reactions based on initial values measured at 1 $m$ meter depth intervals for Upper Mystic Lake in Winchester, MA. The model uses theoretical energetic yields based on redox potentials to determine the order in which nutrients become depleted. \n",
    "\n",
    "The Mystic Lakes are located at 42.4317° N, 71.1483° W, within the boundaries Winchester, Arlington, and Medford, MA. They both lie 1 $m$ above sea level. The Aberjona River feeds into Upper Mystic Lake, which flows over the Mystic dam into Lower Mystic Lake, which drains out through the Mystic River. Mystic River splits within 100m of the lake and a portion goes off into Alewife Creek, while the rest flows out into Boston Harbor.\n",
    "\n",
    "\n",
    "\n",
    "Both Mystic Lakes are dimictic, and therefore both show strong seasonal stratification that progresses throughout the year. In nature, it is assumed that the most of microbial metabolism is fueled by reduction of molecular oxygen for the oxidation of organic carbon. The source of this organic carbon is \"particulate biomass, humics from rainwater runoff, and photosynthesis\" <a name=\"ref-3\"/>[(Preheim et al, 2016)](#cite-Preheim2016). The source of the oxygen is also photosynthesis, as well as diffusion out of the atmosphere, both of which happen just beneath the surface in the hypolimnion. As these two chemical species strongly influence the magnitude and rate of microbial metabolism throughout the water column, I would like to investigate the meterological varaibles which modify their distribution.\n",
    "\n",
    "At present, it uses an empirically accurate methods for modelling vertical transport and mixing. However, the General Lake Model <a name=\"ref-4\"/>[(Read et al, 2014)](#cite-Read2014142) is another similar tool which is better grounded in theory and incorporates meteorological data into determinations of vertical transport. The GLM user manual describes itself as \"ideally suited to long-term investigations ranging from months to decades, and for coupling with biogeochemical models to explore the role that stratification and vertical mixing play on lake ecosystem dynamics.\" It is my ultimate objective to perform a sensitivity analysis on the GLM to determine the relative contributions wind, sun, cloud cover, and precipitation on biogeochemistry within the lake."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Objective\n",
    "\n",
    "The independent variables I would like to evaluate include air temperature, inflow salinity, runoff volume, solar radiation, wind fetch, precipitation, and cloud cover. Each of these variables contribute in varying proportions to the degree of heat flux throughout the water column, and thus control vertical mixing. I have a poor understanding of the relative contributions of each factor, as well as the secondary effects that they control. For example, the duration of sunlight will contribute heat directly to the lake surface, however it will also increase phytoplankton-based turbidity and promote evaporation, which have the opposite effect. The appropriate way of dealing with this is to perform a sort of factor screening experiment <a name=\"ref-5\"/>[(Montgomery, 2000)](#cite-Montgomery2000), in which variables are modified either one at a time, and in combinations, and their individual or interaction effects are observed. The inputs to and the outputs from the model are time series of measurements. Artifical time series will be created such that a baseline condition can be established. The baseline or control conditions will be devised such that they most closely resemble the conditions at Mystic Lake at the appropriate time of year. Data from nearby streams are substituted based on the resemblance of any other available data common between the two locations. A minimum run time of 1 year for each instance of the model. \n",
    "\n",
    "Once a baseline is established, real datasets will be introduced for each varaible. Simulation duration will be limited to minimum of two years. This would provide me with a good intuition into the relative contributions of each input into lake mixing. By calculating the effect of each input dataset on particular output variables below \n",
    "\n",
    "I am interested in seeing how much each variable changes the time series produced for *Lake Number*, *total*, *evaporated*, and *rainfall volumes*, *lake surface area*, *benthic light value*, *average significant wave height*, as well as the *temperatures* at each depth. These 8 parameters were selected from the 27 available options because most are likely indicators of a mixing and a subsequent compositional change within the lake. In particular, significant changes in Lake Number during the Winter season would indicate that deep mixing is occuring and dissolved oxygen is being delivered to the lake bed <a name=\"ref-6\"/>[(Robertson & Imberger, 1994)](#cite-Jorg1994). I would also like to know how much volume runoff contributes into a lake, and the required minimum concentration of water-quality nutrients necessary to significantly modify the fluxes of species considered in the model. If either parameter had a significant impact, it would indicate that integration of biogeochemical models with physical & meteorological models is a necessary marriage in order to adequately characterize microbial metabolism in dimictic lakes (at least)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Software Operation\n",
    "\n",
    "Based on testing performed with the example files provided in the `examples_2.2/warmlake/fabm` directory, the minimum required files for operation include:\n",
    "\n",
    "glm2.nml\n",
    "inflow_1_fabm.csv\n",
    "inflow_2_fabm.csv\n",
    "met_daily.csv\n",
    "outflow.csv\n",
    "\n",
    "The outputs created include:\n",
    "\n",
    "\n",
    "outlet_00.csv\n",
    "overflow.csv\n",
    "lake.csv\n",
    "\n",
    "### output.nc \n",
    "\n",
    "This is a classic netCDF3 archive that contains as many dimensions as necessary to run. The following script shows how to extract, unmask, slice, and plot the underlying information. The example data is provided by the GLM's authors. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data Sets\n",
    "\n",
    "*Figure 1*: Depth contours shown here were the source for the storage, elevation, and area relationships. The same map was used by Fricker & Nepf [](cite-#JGRC8090) in their study on seiches in Upper Mystic Lake.\n",
    "\n",
    "<img src=\"morphometry/UpperMysticUnedited.png\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The areas for the 9 depths provided were calculated using ImageJ and scaled based on satellite image data from Google Maps. obtained to define basin area. These curves provide \"a relatively coarse resolution relative to the simulated layers\" and so the model uses interpolation of volumes of thinner layers (0.1 meters) assuming an underlying conical shape. The bottom layer is calculated first according to the formula $ V_{bottom} = A_{bottom} \\frac{h}{3} $ and the remaining points 8 segments are derived from:\n",
    "\n",
    "$$ V_{upper} = V_{lower} + [A_{lower} + 0.5(A_{upper} - A_{lower})] (ℎ_{upper} − ℎ_{lower}) $$\n",
    "\n",
    "Bathymetry.py: used to plot, scale, and interpolate area-depth relationships"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The daily values for various meteorological & river flow are also required in the GLM. The following table shows the data type, frequency, name, frequency, source, and file name in the model.\n",
    "\n",
    "|Name|Source|Input File Name|Units Provided|Units Req'd|Coverage|\n",
    "|---|---|---|---|---|---|\n",
    "|Inflow Volume/Discharge|USGS|`inflows.csv`|${ft}^3/second$|$ml/day$|10/1/2011-3/19/2016|\n",
    "|Outflow Volume/Discharge|USGS|`outflows.csv`|${ft}^3/second$|$ml/day$|9/23/2015-3/19/2016|\n",
    "|Inflow Salinity Spec. conductance|USGS|`inflows.csv`|$\\frac{\\mu S}{cm}$ @ 25°C|$mg/L$|11/15/2011-3/19/2016|\n",
    "|Inflow Temperature|USGS|`inflows.csv`|°C|°C|11/15/2011-3/19/2016|\n",
    "|Shortwave Radiation/Irradiance|CERES|`met.csv`|$W/m^2$|$W/m^2$|1/1/2010-1/1/2014|\n",
    "|Longwave Radiation/Irradiance|CERES|`met.csv`|$W/m^2$|$W/m^2$|1/1/2010-1/1/2014|\n",
    "|Cloud Cover (Fraction)|LAADS|`met.csv`|?|Fraction|1/1/2010-1/1/2014|\n",
    "|Air Temperature|CERES|`met.csv`|K|°C|1/1/2010-1/1/2014|\n",
    "|Relative Humidity|CERES|`met.csv`|%|%|1/1/2010-1/1/2014|\n",
    "|Wind Speed|GHCN|`met.csv`|24 hr|$m/s$|1/1/2000-1/1/2015|\n",
    "|Rainfall|GHCN|`met.csv`|$\\frac{mm}{10}$|$m/day$|1/1/2000-1/1/2015|\n",
    "|Snowfall|GHCN|`met.csv`|$\\frac{mm}{10}$|$m/day$|1/1/2000-1/1/2015|"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Weather Data\n",
    "\n",
    "### Global Historical Climatology Network Data\n",
    "\n",
    "Some meterological data was sourced from a [GHCN](https://gis.ncdc.noaa.gov/geoportal/catalog/search/resource/details.page?id=gov.noaa.ncdc:C00861) station 8.35 miles East South East at Boston Logan Airport. The GCHN data was the easiest to read in and contained only 2 missing values across all of 2000-2015. It contains daily values on **wind speed, temperature min/max, precipitation, & snowfall.** It is called as an object with an associated Pandas dataframe.\n",
    "\n",
    "```python\n",
    "met_fn = \"BostonLoganAirportWeather.csv\"\n",
    "met_path = os.path.join(mainFolder, 'weatherData', met_fn)\n",
    "BOS_weather = LakeModel.GHCN_weather_data(met_path)\n",
    "# Makes the column names more readable\n",
    "BOS_weather.clean_columns()\n",
    "```\n",
    "TODO: plot that data here\n",
    "\n",
    "### Clouds and the Earth's Radiant Energy System Data\n",
    "\n",
    "[NASA Ceres](http://ceres.larc.nasa.gov/) data is derived from NASAs Earth Observing System, MODIS, and VIRS projects. Two data products were utilized from this database. The first was the Level 3 SYN-deg product which provided daily data at a spatial resolution of 20 sq. km for two location within the requested bounding box (42.5 N, 71.5 W & 70.5 W). The actual data types contained include **short & long wave radiation** measured in both directions provided at the surface, the **air temperature** measured at the skin surface (?), and **cloud cover**. Access protocol is shown here:\n",
    "\n",
    "```python\n",
    "ceres_fn2 = \"CERES_SYN1deg-Day_200508-201406.nc\"\n",
    "ceres_SYN1 = os.path.join(mainFolder, 'weatherData', ceres_fn2)\n",
    "CERES_data = LakeModel.CERES_nc(ceres_SYN1, \"4\")\n",
    "\n",
    "# df1 contains more accurate data for 42.5 N, 71.5 W. df2 applies to 42.5 N, 70.5 W. \n",
    "for i in CERES_data.df1.columns:\n",
    "    print i\n",
    "\n",
    "    sfc_comp_sw-up_all_daily\n",
    "    sfc_comp_sw-down_all_daily\n",
    "    sfc_comp_lw-up_all_daily\n",
    "    sfc_comp_lw-down_all_daily\n",
    "    cldarea_total_daily\n",
    "    aux_skint_daily\n",
    "    aux_surfpress_daily\n",
    "    \n",
    "```\n",
    "\n",
    "The cloud area data was compared to the less mature Level 2 SSF product from CERES and the OMMYDCLD product obtained from the [Goddard Earth Sciences Data and Information Service Center](http://daac.gsfc.nasa.gov/), which is based mostly on the MODIS06_L2 product. The cloud cover data within the SYN-deg product better resembled the OMMYDCLD data (R = 0.71, R=0.73), compared to the less processed Level SSF (R = 0.15, 0.21), despite the fact that it is the same information, subject to spatial and temportal smoothing. A few months of the SSF & Giovanni signals are shown below to demonstrate the phase agreement.\n",
    "\n",
    "<img src=\"donePlots/CERES_SSFvOMMYCLD0.png\", width=800, height=600>\n",
    "\n",
    "The SSF data object is still required for use of the humidity, although other metrics are contained. Access and preprocessing are shown here:\n",
    "\n",
    "```python\n",
    "ceres_fn = \"CERES_SSF_XTRK-MODIS_Edition3A_Subset_2010010102-2014010116.nc\"\n",
    "ceres_ssf = os.path.join(mainFolder, 'weatherData', ceres_fn)\n",
    "CERES_SSF = LakeModel.CERES_nc(ceres_ssf, \"3\")\n",
    "Geo_D = CERES_SSF.ceres_df.resample('D').mean()\n",
    "print Geo_D.columns\n",
    "    Index([u'temp', u'windU', u'windV', u'humidity', u'cloud_frac', u'LW_rad',\n",
    "           u'SW_rad'], dtype='object')\n",
    "```\n",
    "\n",
    "The SYN data spans 2005-01-01 to 2014-06-30 and the SSF data covers 2010-01-01 to 2014-01-01. It is interesting to note that while General Lake Model only uses the daily average wind speed, there are other more informative measures that we obtained incidentally. In addition to the average daily wind speed (`AWND`), the fastest 2-minute/5-second wind speed & direction (`WSF2` and `WSF5`) are also present. The CERES data set goes even further in that it provides surface Wind U-vector and surface wind V-vector. As these measures make explicit the direction of the force of wind, it is possible to model the effect Mystic Lake's asymetrical morphometry has on sieche propogation.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Water Data\n",
    "\n",
    "Data about inflows and outflows was obtained from the USGS Water Data site, using a search query for sites with >100 observations, within 42.50 to 42.35 N and -71.35 to -71.00 W. No public data is available for the Mystic and Aberjona rivers on water temperature & conductance data. Such data *does* exist for the Hobbs Brook & Fresh Pond Reservoirs, both less than 11.06 km away and within the same watershed. Although not used by this work, there is also bathymetry for Hobbs Brook, a record of Groundwater data very near Mystic Lake, and monthly values for water quality maintained by the [Mystic River Watershed Assoc.](http://mysticriver.org/baseline/).\n",
    "\n",
    "The `waterdata/` folder contains data from the following sites: \n",
    "\n",
    " - [Aberonja River](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01102500)\n",
    "     - Period: 1939-04-16 to 2016-03-31\n",
    "     - Data Types: Discharge \n",
    "     - Description: Main inflowing river\n",
    " - [Mystic River](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01103010)\n",
    "     - Period: 2015-09-23 to 2016-06-30\n",
    "     - Data Types: Discharge, Water Temperature (15 min interval) \n",
    "     - Description: Note, there are two files for this data because of the different time interval \n",
    " - [Mystic River Rt 16 @ Medford, MA](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01103040)\n",
    "     - Period: 2016-03-01 to 2016-06-30\n",
    "     - Data Types: Discharge\n",
    "     - Description: This is another site farther downstream     \n",
    " - [Alewife Brook near Arlington, MA](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01103025)\n",
    "     - Period: 2005-10-01 to 2016-06-29\n",
    "     - Data Types: Discharge\n",
    "     - Description: Lat 42o24'25\", long 71o08'04\". This feeds into the Mystic River quite close to where Mystic Lake does the same \n",
    " - [Fresh Pond](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=422302071083801)\n",
    "     - Period: 2005-10-02 to 2016-06-29\n",
    "     - Data Types: Temperature (Air), Temperature (Water), Specific conductance, Precipitation\n",
    "     - Description: Very nearby, but smaller reservoir \n",
    " - [Hobbs Brook Source @ Mill St.](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01104405)\n",
    "    - One of four tributaries that feeds Cambridge Reservoir\n",
    "    - Period: 2011-11-15 to 2015-09-06 \n",
    "    - Discharge, Conductance, and Water Temperature\n",
    " - [Cambridge Res., Unnamed Tributary 1](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01104410)\n",
    "     - Period: 2011-11-16 to 2015-08-07\n",
    "     - Data Types: Discharge, Conductance, and Water Temperature\n",
    "     - Description: One of three tributaries that feeds Cambridge Reservoir\n",
    " - [Cambridge Res., Unnamed Tributary 2](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01104415)\n",
    "     - Period: 2004-10-01 to 2016-06-29\n",
    "     - Data Types: Discharge, Conductance, and Water Temperature\n",
    "     - Description: One of three tributaries that feeds Cambridge Reservoir\n",
    " - [Cambridge Res., Unnamed Tributary 3](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01104420)\n",
    "     - Period: 2012-03-26 to 2016-06-29\n",
    "     - Data Types: Discharge, Conductance, and Water Temperature\n",
    "     - Description: One of three tributaries that feeds Cambridge Reservoir\n",
    " - [Hobbs Brook below Cambridge](http://waterdata.usgs.gov/nwis/inventory?agency_code=USGS&site_no=01104430)\n",
    "      - 2006-10-01 to 2015-10-07 \n",
    "      - Period: Discharge, Conductance, Air Temp, Water Temp, Precipitation\n",
    "\n",
    "There is a good map and a description of the directional flow relationships between some of the water bodies mentioned below at [the Cambridge municipality website](https://www.cambridgema.gov/Water/watershedmanagementdivision/wheredoesourdrinkingwatercomefrom). To determine whether some of these sources could be used as in place of data specific to the site, the following figures were created. The first is a comparison between Hobbs Brook below Cambridge Res. and the Fresh Pond Res. The strong disparity in seasonal conductance may be attributable to the recharge effect described at the previous link. \n",
    "\n",
    "<img src=\"donePlots/HobbsandFresh.png\", width=800, height=600>\n",
    "\n",
    "The second plot is seasonally-aligned values for surface water temperatures at Fresh Pond and Hobbs Brook in 2012, compared to those measured by us at Mystic Lake in 2007.\n",
    "\n",
    "<img src=\"donePlots/MysticLakeHobbsandFreshTemperatures.png\", width=800, height=600>\n",
    "\n",
    "- (CERES Wind speed, GHCN Wind Speed)\n",
    "- (Aberjona Precipitation, GHCN Precipitation)\n",
    "\n",
    "** Milestone 4**: Optimize Lake Specific Model Parameters<br>\n",
    "The following initial conditions may also be modified as they are \"lake specific\": *Stream-bed drag*, *seepage rate*, *critical area*, and the *extinction coefficient* for PAR radiation. The criticial area is \"critical area below which wind sheltering may occur.\" I am not sure how to calculate appropriate values for these, so testing the effect of deviations from their defaults may prove interesting.\n",
    "\n",
    "** Milestone 5**: Run model with baseline and experimental conditions <br>\n",
    "Write script to automate iterative evaluations of `GLM` with each combination of input conditions, record dependent variables, and, once complete, test results for statistical significance. \n",
    "\n",
    "** Milestone 6**: Modify MATLAB implementation (shown below) to handle any significant meterological forcing inputs <br>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The relationship between discharge of the Aberjona and the Mystic Rivers for the short overlapping time period is not entirely clear, as shown by the next figure. There is clearly a relationship, but it appears to be highly buffered by the volume of the lake and it is not clear to me whether there is enough data here to model that effect. \n",
    "\n",
    "<img src=\"donePlots/DammedDischarge.png\", width=800, height=600>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mystic_d.df.interpolate(inplace=True)\n",
    "inflow_s, mystic_ds = LakeModel.subsectbydate_2(inflow.df, mystic_d.df)\n",
    "ground_s, _ = LakeModel.subsectbydate_2(GroundW.df, mystic_d.df)\n",
    "plt.plot(inflow_s.index, inflow_s.values, label='inflow')\n",
    "plt.plot(mystic_ds.index, mystic_ds.values, label='outflow')\n",
    "plt.plot(ground_s.index, ground_s.values, label='groundwater')\n",
    "plt.legend(loc='best')\n",
    "if not os.path.exists(\"../donePlots/Mystic&Aberjona_Discharge.png\"):\n",
    "    plt.savefig(\"../donePlots/Mystic&Aberjona_Discharge.png\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- (CERES Wind speed, GHCN Wind Speed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'CERES_data' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-10-109b472e84f7>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m#temperature comparison\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mt1\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mCERES_data\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdf1\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0maux_skint_daily\u001b[0m \u001b[0;31m#2005-2014\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m \u001b[0mt2\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mBOS_weather\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdf\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mt_max\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mt3\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mBOS_weather\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdf\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mt_min\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mt4\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mHobbs_o\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdf\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'Temperature, air, degrees Celsius (Mean)'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdropna\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'CERES_data' is not defined"
     ]
    }
   ],
   "source": [
    "#temperature comparison\n",
    "t1 = CERES_data.df1.aux_skint_daily #2005-2014\n",
    "t2 = BOS_weather.df.t_max\n",
    "t3 = BOS_weather.df.t_min\n",
    "t4 = Hobbs_o.df['Temperature, air, degrees Celsius (Mean)'].dropna()\n",
    "t5 = Hobbs_o.df['Temperature, air, degrees Celsius (Minimum)'].dropna()\n",
    "plt.plot(t1.index, t1.values, label='Skin Temperature (CERES)')\n",
    "plt.plot(t2.index, t2.values, label='Max Temp (recorded at BOS airport)')\n",
    "plt.plot(t3.index, t3.values, label='Min Temp (recorded at BOS airport)')\n",
    "plt.plot(t4.index, t4.values, label='Mean Temp (recorded at Hobbs Brook)')\n",
    "plt.plot(t5.index, t5.values, label='Min Temp (recorded at Hobbs Brook)')\n",
    "plt.xlim([t4.index[1188], t4.index[1188+365]])\n",
    "plt.legend(loc='best')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Current Implementation \n",
    "\n",
    "The Mystic Lake model solves an equation called `flux` using the `ODE15s` solver. This function is taken as an argument by the ODE solver, which feeds this function the concentration matrix (flattened into a concentration vector). The function computes the rates of reactions, diffusions, and precipitations and outputs the fluxes for each metabolite at each depth."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "\n",
    "```matlab\n",
    "% twiddle because this is time-independent\n",
    "function [conc_fluxes] = flux(~, concs_vector)\n",
    "\n",
    "    % These are the matrices containing the concentrations and fluxes of each species\n",
    "    % There are 10 species (columns) and 17 layers (1 meter depth sections)\n",
    "    concs = reshape(concs_vector, [n_x, n_species]);\n",
    "    conc_fluxes = zeros(n_x, n_species);\n",
    "\n",
    "    %... some code is omitted ...\n",
    "    \n",
    "    % loop from surface to bottom\n",
    "    for x = 1: n_x\n",
    "        %% ... skip down to the relevant part ....\n",
    "\n",
    "        % diffusion is based on a D_plus and a D_minus coefficients for all species\n",
    "        % The following lines state that: ...\n",
    "        if x > 1\n",
    "            conc_fluxes(x, :) = conc_fluxes(x, :) + D_plus .* concs(x - 1, :) - D_minus .* concs(x, :);\n",
    "            % this line adds what diffuses down from above and subtracts what diffuses up from this layer\n",
    "        end\n",
    "\n",
    "        if x < n_x\n",
    "            conc_fluxes(x, :) = conc_fluxes(x, :) - D_plus .* concs(x, :) + D_minus .* concs(x + 1, :);\n",
    "           \t% this line subtracts what diffuses down from this layer, and addes what diffuses up from below\n",
    "        end\n",
    "    end\n",
    "    % .... more code is omitted ...\n",
    "end\n",
    "\n",
    "% Call to solver function & options provided to it.\n",
    "options = odeset('NonNegative', 1: n_total);\n",
    "[time_slices, y] = ode15s(@flux, linspace(0.0, t_max, n_time_slices), concs0_vector, options);\n",
    "```\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "<!--bibtex\n",
    "\n",
    "@article{Hunter199853,\n",
    "title = \"Kinetic modeling of microbially-driven redox chemistry of subsurface environments: coupling transport, microbial metabolism and geochemistry \",\n",
    "journal = \"Journal of Hydrology \",\n",
    "volume = \"209\",\n",
    "number = \"1–4\",\n",
    "pages = \"53 - 80\",\n",
    "year = \"1998\",\n",
    "note = \"\",\n",
    "issn = \"0022-1694\",\n",
    "doi = \"http://dx.doi.org/10.1016/S0022-1694(98)00157-7\",\n",
    "url = \"http://www.sciencedirect.com/science/article/pii/S0022169498001577\",\n",
    "author = \"Kimberley S. Hunter and Yifeng Wang and Philippe Van Cappellen\",\n",
    "keywords = \"Biodegradation kinetics\",\n",
    "keywords = \"Redox chemistry\",\n",
    "keywords = \"Reactive transport\",\n",
    "keywords = \"Microbial activity\",\n",
    "keywords = \"Subsurface environment \"\n",
    "}\n",
    "@article{Preheim2016,\n",
    "title = \"Surveys, simulation, and single-cell assays link function to phylogeny in a natural ecosystem\",\n",
    "journal = \"In Review\",\n",
    "year = \"2016\",\n",
    "author = \"Sarah P. Preheim, Scott W. Oleson, Sarah J. Spencer, Arne Materna, Charuleka Varadharajan, Matthew Blackburn, Jonatham Friedman, Jorge Rodriguiez, Harold Hemond, and Eric J. Alm\"\n",
    "}\n",
    "@article{Read2014142,\n",
    "title = \"Simulating 2368 temperate lakes reveals weak coherence in stratification phenology \",\n",
    "journal = \"Ecological Modelling \",\n",
    "volume = \"291\",\n",
    "number = \"\",\n",
    "pages = \"142 - 150\",\n",
    "year = \"2014\",\n",
    "note = \"\",\n",
    "issn = \"0304-3800\",\n",
    "doi = \"http://dx.doi.org/10.1016/j.ecolmodel.2014.07.029\",\n",
    "url = \"http://www.sciencedirect.com/science/article/pii/S0304380014003664\",\n",
    "author = \"Jordan S. Read and Luke A. Winslow and Gretchen J.A. Hansen and Jamon Van Den Hoek and Paul C. Hanson and Louise C. Bruce and Corey D. Markfort\",\n",
    "keywords = \"Limnology\",\n",
    "keywords = \"Climate change\",\n",
    "keywords = \"Water temperature\",\n",
    "keywords = \"Lake temperature modeling\",\n",
    "keywords = \"Coherence\",\n",
    "keywords = \"Stratification \"\n",
    "}\n",
    "@article{Montgomery2000,\n",
    "title = \"Design and Analysis of Experiments. 7th edition\",\n",
    "journal = \"John Wiley & Sons, Inc\",\n",
    "author = \"Douglas C. Montgomery\",\n",
    "year = \"2000\",\n",
    "url = \"http://www.wiley.com/WileyCDA/WileyTitle/productCd-EHEP002024.html\"\n",
    "}\n",
    "@article{JGRC8090,\n",
    "author = {Fricker, P. D. and Nepf, H. M.},\n",
    "title = {Bathymetry, stratification, and internal seiche structure},\n",
    "journal = {Journal of Geophysical Research: Oceans},\n",
    "volume = {105},\n",
    "number = {C6},\n",
    "issn = {2156-2202},\n",
    "url = {http://dx.doi.org/10.1029/2000JC900060},\n",
    "doi = {10.1029/2000JC900060},\n",
    "pages = {14237--14251},\n",
    "year = {2000},\n",
    "}\n",
    "@article {Jorg1994,\n",
    "author = {Robertson, Dale M. and Imberger, Jörg},\n",
    "title = {Lake Number, a Quantitative Indicator of Mixing Used to Estimate Changes in Dissolved Oxygen},\n",
    "journal = {Internationale Revue der gesamten Hydrobiologie und Hydrographie},\n",
    "volume = {79},\n",
    "number = {2},\n",
    "publisher = {Akademie Verlag, Berlin},\n",
    "issn = {1522-2632},\n",
    "url = {http://dx.doi.org/10.1002/iroh.19940790202},\n",
    "doi = {10.1002/iroh.19940790202},\n",
    "pages = {159--176},\n",
    "keywords = {Lakes, mixing, dissolved oxygen},\n",
    "year = {1994},\n",
    "}\n",
    "-->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# References\n",
    "\n",
    "<a name=\"cite-Hunter199853\"/>1. Kimberley S. Hunter and Yifeng Wang and Philippe Van Cappellen. 1998. **Kinetic modeling of microbially-driven redox chemistry of subsurface environments: coupling transport, microbial metabolism and geochemistry **. _Journal of Hydrology _. [URL](http://www.sciencedirect.com/science/article/pii/S0022169498001577)\n",
    "\n",
    "<a name=\"cite-Preheim2016\"/>2. Sarah P. Preheim, Scott W. Oleson, Sarah J. Spencer, Arne Materna, Charuleka Varadharajan, Matthew Blackburn, Jonatham Friedman, Jorge Rodriguiez, Harold Hemond, and Eric J. Alm. 2016. **Surveys, simulation, and single-cell assays link function to phylogeny in a natural ecosystem**. _In Review_.\n",
    "\n",
    "<a name=\"cite-Read2014142\"/>3. Jordan S. Read and Luke A. Winslow and Gretchen J.A. Hansen and Jamon Van Den Hoek and Paul C. Hanson and Louise C. Bruce and Corey D. Markfort. 2014. **Simulating 2368 temperate lakes reveals weak coherence in stratification phenology **. _Ecological Modelling _. [URL](http://www.sciencedirect.com/science/article/pii/S0304380014003664)\n",
    "\n",
    "<a name=\"cite-Montgomery2000\"/>4. Douglas C. Montgomery. 2000. **Design and Analysis of Experiments. 7th edition**. _John Wiley & Sons, Inc_. [URL](http://www.wiley.com/WileyCDA/WileyTitle/productCd-EHEP002024.html)\n",
    "\n",
    "<a name=\"cite-Jorg1994\"/>5. Robertson, Dale M. and Imberger, Jörg. 1994. **Lake Number, a Quantitative Indicator of Mixing Used to Estimate Changes in Dissolved Oxygen**. _Internationale Revue der gesamten Hydrobiologie und Hydrographie_. [URL](http://dx.doi.org/10.1002/iroh.19940790202)\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  },
  "latex_metadata": {
   "affiliation": "Johns Hopkins University",
   "author": "Keith Arora-Williams"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}