RRS Final Concatination

RRS Final Concatination#

Now that SeaBASS and GLORIA RRS data has been organized and standardized into the same format, they can be combined into a single dataset

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import cartopy
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import geopandas as gpd
from datetime import datetime
import os
from matplotlib import ticker
import datetime as dt
import plotly.express as px
import cmocean as cm
import cmocean.cm as cmo
import matplotlib.gridspec as gridspec
import time
import matplotlib.ticker as mticker
gloria = pd.read_excel(r'C:\Users\gianna.milton\Documents\Python\Coastal_chl_final\GLORIA_rrs_na.xlsx')
gloria['datetime'] = pd.to_datetime(gloria['datetime'])
gloria['source']='GLORIA'
gloria = gloria[gloria['datetime'] >= '2000-01-01'] #only want data from 2000 on for this algorithm 
---------------------------------------------------------------------------
KeyboardInterrupt                         Traceback (most recent call last)
Cell In[2], line 1
----> 1 gloria = pd.read_excel(r'C:\Users\gianna.milton\Documents\Python\Coastal_chl_final\GLORIA_rrs_na.xlsx')
      2 gloria['datetime'] = pd.to_datetime(gloria['datetime'])
      3 gloria['source']='GLORIA'

File ~\AppData\Local\anaconda3\Lib\site-packages\pandas\io\excel\_base.py:508, in read_excel(io, sheet_name, header, names, index_col, usecols, dtype, engine, converters, true_values, false_values, skiprows, nrows, na_values, keep_default_na, na_filter, verbose, parse_dates, date_parser, date_format, thousands, decimal, comment, skipfooter, storage_options, dtype_backend, engine_kwargs)
    502     raise ValueError(
    503         "Engine should not be specified when passing "
    504         "an ExcelFile - ExcelFile already has the engine set"
    505     )
    507 try:
--> 508     data = io.parse(
    509         sheet_name=sheet_name,
    510         header=header,
    511         names=names,
    512         index_col=index_col,
    513         usecols=usecols,
    514         dtype=dtype,
    515         converters=converters,
    516         true_values=true_values,
    517         false_values=false_values,
    518         skiprows=skiprows,
    519         nrows=nrows,
    520         na_values=na_values,
    521         keep_default_na=keep_default_na,
    522         na_filter=na_filter,
    523         verbose=verbose,
    524         parse_dates=parse_dates,
    525         date_parser=date_parser,
    526         date_format=date_format,
    527         thousands=thousands,
    528         decimal=decimal,
    529         comment=comment,
    530         skipfooter=skipfooter,
    531         dtype_backend=dtype_backend,
    532     )
    533 finally:
    534     # make sure to close opened file handles
    535     if should_close:

File ~\AppData\Local\anaconda3\Lib\site-packages\pandas\io\excel\_base.py:1616, in ExcelFile.parse(self, sheet_name, header, names, index_col, usecols, converters, true_values, false_values, skiprows, nrows, na_values, parse_dates, date_parser, date_format, thousands, comment, skipfooter, dtype_backend, **kwds)
   1576 def parse(
   1577     self,
   1578     sheet_name: str | int | list[int] | list[str] | None = 0,
   (...)
   1596     **kwds,
   1597 ) -> DataFrame | dict[str, DataFrame] | dict[int, DataFrame]:
   1598     """
   1599     Parse specified sheet(s) into a DataFrame.
   1600 
   (...)
   1614     >>> file.parse()  # doctest: +SKIP
   1615     """
-> 1616     return self._reader.parse(
   1617         sheet_name=sheet_name,
   1618         header=header,
   1619         names=names,
   1620         index_col=index_col,
   1621         usecols=usecols,
   1622         converters=converters,
   1623         true_values=true_values,
   1624         false_values=false_values,
   1625         skiprows=skiprows,
   1626         nrows=nrows,
   1627         na_values=na_values,
   1628         parse_dates=parse_dates,
   1629         date_parser=date_parser,
   1630         date_format=date_format,
   1631         thousands=thousands,
   1632         comment=comment,
   1633         skipfooter=skipfooter,
   1634         dtype_backend=dtype_backend,
   1635         **kwds,
   1636     )

File ~\AppData\Local\anaconda3\Lib\site-packages\pandas\io\excel\_base.py:778, in BaseExcelReader.parse(self, sheet_name, header, names, index_col, usecols, dtype, true_values, false_values, skiprows, nrows, na_values, verbose, parse_dates, date_parser, date_format, thousands, decimal, comment, skipfooter, dtype_backend, **kwds)
    775     sheet = self.get_sheet_by_index(asheetname)
    777 file_rows_needed = self._calc_rows(header, index_col, skiprows, nrows)
--> 778 data = self.get_sheet_data(sheet, file_rows_needed)
    779 if hasattr(sheet, "close"):
    780     # pyxlsb opens two TemporaryFiles
    781     sheet.close()

File ~\AppData\Local\anaconda3\Lib\site-packages\pandas\io\excel\_openpyxl.py:615, in OpenpyxlReader.get_sheet_data(self, sheet, file_rows_needed)
    613 data: list[list[Scalar]] = []
    614 last_row_with_data = -1
--> 615 for row_number, row in enumerate(sheet.rows):
    616     converted_row = [self._convert_cell(cell) for cell in row]
    617     while converted_row and converted_row[-1] == "":
    618         # trim trailing empty elements

File ~\AppData\Local\anaconda3\Lib\site-packages\openpyxl\worksheet\_read_only.py:85, in ReadOnlyWorksheet._cells_by_row(self, min_col, min_row, max_col, max_row, values_only)
     77 with self._get_source() as src:
     78     parser = WorkSheetParser(src,
     79                              self._shared_strings,
     80                              data_only=self.parent.data_only,
     81                              epoch=self.parent.epoch,
     82                              date_formats=self.parent._date_formats,
     83                              timedelta_formats=self.parent._timedelta_formats)
---> 85     for idx, row in parser.parse():
     86         if max_row is not None and idx > max_row:
     87             break

File ~\AppData\Local\anaconda3\Lib\site-packages\openpyxl\worksheet\_reader.py:156, in parse()

File ~\AppData\Local\anaconda3\Lib\xml\etree\ElementTree.py:1238, in iterparse.<locals>.iterator(source)
   1236 yield from pullparser.read_events()
   1237 # load event buffer
-> 1238 data = source.read(16 * 1024)
   1239 if not data:
   1240     break

File ~\AppData\Local\anaconda3\Lib\zipfile\__init__.py:989, in ZipExtFile.read(self, n)
    987 self._offset = 0
    988 while n > 0 and not self._eof:
--> 989     data = self._read1(n)
    990     if n < len(data):
    991         self._readbuffer = data

File ~\AppData\Local\anaconda3\Lib\zipfile\__init__.py:1065, in ZipExtFile._read1(self, n)
   1063 elif self._compress_type == ZIP_DEFLATED:
   1064     n = max(n, self.MIN_READ_SIZE)
-> 1065     data = self._decompressor.decompress(data, n)
   1066     self._eof = (self._decompressor.eof or
   1067                  self._compress_left <= 0 and
   1068                  not self._decompressor.unconsumed_tail)
   1069     if self._eof:

KeyboardInterrupt: 
seabass = pd.read_excel(r'C:\Users\gianna.milton\Documents\Python\Coastal_chl_final\SB_rrs_na.xlsx')
seabass = seabass.rename(columns={'identifier_product_doi':'DOI_url'})
seabass=seabass[['datetime', 'lon', 'lat', 'DOI_url', 'affiliations', 'investigators', 'contact', 'experiment', 'cruise',
       'station', 'rrs', 'wavelength', 'depth']]
seabass['source']='SeaBASS'
seabass['datetime'] = pd.to_datetime(seabass['datetime'])

Since some GLORIA data is duplicated on SeaBASS, remove any duplicates from the SeaBASS dataset for better program cohesion

#remove GLORIA data from Seabass 
gloria_sb = gloria[gloria['DOI_url'] !='https://doi.pangaea.de/10.1594/PANGAEA.948492'] #only gather data that has seabass doi
gloria_sb=gloria_sb.rename(columns={"DOI_url": "doi_url"})
gloria_sb=gloria_sb[['datetime', 'doi_url','lat', 'lon','wavelength','rrs']]

seabass_test=seabass.copy()
seabass_test['doi_url'] = seabass_test['DOI_url'].str.lower()

test2 = pd.merge(gloria_sb, seabass_test, on=['datetime', 'doi_url','lat', 'lon','wavelength','rrs'], how='right',indicator=True).reset_index(drop=True)
#seabass data has some wavelengths that gloria does not i.e., 10.5067/seabass/2009oct_chesapeake/data001 has wavelengths above 800 on seabass not on gloria
seabass =  test2[test2['_merge'] != 'both'] #remove rows where _merge has both
seabass=seabass[['datetime', 'lat', 'lon', 'wavelength', 'rrs','source', 'DOI_url', 'affiliations', 'investigators', 'contact', 'experiment', 
                 'cruise','station', 'depth']]
dfs=[seabass,gloria]
all_rrs = pd.concat(dfs).reset_index(drop=True)
all_vars = ['source', 'datetime', 'lon', 'lat', 'depth', 'rrs', 'wavelength',  'experiment', 'DOI_url', 'affiliations', 'investigators', 'contact', 
            'cruise', 'station']
all_rrs = all_rrs.drop_duplicates(subset=all_vars, keep='first')
all_rrs['datetime'] = pd.to_datetime(all_rrs['datetime'])

Currently, the combined dataset is in long format, where each unique wavelength is in a seperate column from the rrs value. This means that samples have many repeated metadata rows to account for multiple wavelengths. To make the ID column for each unique sample, the dataframe has to be transformed into ‘wide’ format to properly count duplicates and create unique IDs.

def long_to_wide(df):
    """
    Transforms long format to wide format (rrs_###) without averaging duplicates,
    preserving rows even if they have missing metadata (NaNs).
    """
    id_vars = ['source', 'datetime', 'lon', 'lat', 'depth', 'experiment', 'DOI_url', 'affiliations', 'investigators', 'contact', 'cruise', 'station']
    # temporarily fill NaNs with a string so pivot_table doesn't drop them
    df[id_vars] = df[id_vars].fillna('MISSING_DATA')
    
    #group by the full id_vars list to ensure the counter perfectly aligns
    df['temp_counter'] = df.groupby(id_vars + ['wavelength']).cumcount()
    df_wide = df.pivot_table(index=id_vars + ['temp_counter'], columns='wavelength',  values='rrs')
    new_column_names = [f"rrs_{str(col)}" for col in df_wide.columns]
    df_wide.columns = new_column_names
    df_wide = df_wide.reset_index()
    df_wide = df_wide.drop(columns=['temp_counter'])
    return df_wide

test = long_to_wide(all_rrs)
test['temp_exp'] = test['experiment'].str.replace('_', '', regex=False)
test['temp_exp'] = test['temp_exp'].str.replace('-', '', regex=False)
test['temp_exp'] = test['temp_exp'].str.replace(' ', '', regex=False)
test['temp_exp'] = test['temp_exp'].str.replace('(', '', regex=False)
test['temp_exp'] = test['temp_exp'].str.replace(')', '', regex=False)
test['temp_exp'] = test['temp_exp'].str.replace('[', '', regex=False)
test['temp_exp'] = test['temp_exp'].str.replace(']', '', regex=False)

test['ID_code'] = (test['source'].astype(str) + '_' + test['temp_exp'].astype(str) + '_' + test['datetime'].dt.strftime('%Y%m%d-%H%M%S').astype(str) + '_' + 
    test['lat'].astype(str) + '_' + test['lon'].astype(str) + '_' + test['depth'].astype(str) + 'm')
test['ID_code'] = test['ID_code'] + '_' + test.groupby('ID_code').cumcount().astype(str)
#turn back into rrs column adn wavelength column
rrs_cols = [col for col in test.columns if col.startswith('rrs_')]
id_vars_melt = ['source', 'datetime', 'lon', 'lat', 'depth', 'experiment', 'DOI_url', 'affiliations', 'investigators', 'contact', 'cruise', 'station', 'ID_code']

df_long = test.melt(id_vars=id_vars_melt, value_vars=rrs_cols,var_name='raw_wavelength',value_name='rrs')
df_long['raw_wavelength'] = df_long['raw_wavelength'].str.replace('rrs_', '')
df_long['wavelength'] = pd.to_numeric(df_long['raw_wavelength'])
df_long = df_long.drop(columns=['raw_wavelength'])
df_long = df_long.dropna(subset=['rrs'])
df_long[id_vars_melt] = df_long[id_vars_melt].replace('MISSING_DATA', np.nan)
df_long=df_long.reset_index(drop=True)
#df_long.to_excel('all_rrs.xlsx', index = False)

Now all RRS data is standardized and combined into one dataset.

Plots#

rrs = pd.read_excel(r'C:\Users\gianna.milton\Documents\Python\Coastal_chl_final\all_rrs.xlsx')
#first, turn rrs data long format so that there is 1 uinique id for each row / group of wavelengths
def long_to_wide(df):
    """
    Transforms long format to wide format (rrs_###) without averaging duplicates,
    preserving rows even if they have missing metadata (NaNs).
    """
    id_vars = ['source', 'datetime', 'lon', 'lat', 'depth', 'experiment', 'DOI_url','affiliations', 'investigators', 'contact', 'cruise', 'station',
               'ID_code']
    df_temp = df.copy()
    # temporarily fill NaNs with a string so pivot_table doesn't drop them
    df_temp[id_vars] = df_temp[id_vars].fillna('MISSING_DATA')
    
    #group by the full id_vars list to ensure the counter perfectly aligns
    df_temp['temp_counter'] = df_temp.groupby(id_vars + ['wavelength']).cumcount()
    df_wide = df_temp.pivot_table(index=id_vars + ['temp_counter'], columns='wavelength',  values='rrs')
    new_column_names = [f"rrs_{str(col)}" for col in df_wide.columns]
    df_wide.columns = new_column_names
    df_wide = df_wide.reset_index()
    df_wide = df_wide.drop(columns=['temp_counter'])
    return df_wide

rrs_long = long_to_wide(rrs)
rrs_long=rrs_long.replace('MISSING_DATA', np.nan)
category_counts = rrs['source'].value_counts() 
plt.figure(figsize=(5, 5)) # Optional: set the figure size
category_counts.plot.pie(autopct='%1.1f%%', startangle=90, cmap='tab10')
plt.axis('equal') 
plt.show()
_images/32d68140b4da58950966619bc9026bec1da0855bf8c5f629f5f6bf9b6db9f1d5.png
year_test=rrs.copy()
year_test['datetime'] = pd.to_datetime(year_test['datetime'])
year_test['year'] = year_test['datetime'].dt.year
grouped = year_test.groupby(['year', 'source']).size().reset_index(name='DataPoints')

# Create bar chart
fig = px.bar(grouped, x='year', y='DataPoints', color='source', title='Yearly distribution all data',
    labels={'year': 'Year', 'DataPoints': 'Number of Data Points', 'metadata': 'Metadata'},)
fig.update_xaxes(range=[1999,2027])
fig.update_layout(barmode='stack')  # ensures stacking
fig.show()
_images/6a2169b95773db073350c327dde5066b25a0c0e1cce8f2273c16988251d98e1a.png
from matplotlib.colors import LogNorm # Important for high-variance data
fig = plt.figure(figsize=(15, 10))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
ax.add_feature(cfeature.LAND)
ax.add_feature(cfeature.OCEAN)
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.BORDERS)
ax.add_feature(cfeature.STATES)
hb = ax.hexbin(rrs.lon, rrs.lat, gridsize=75, cmap='inferno_r', mincnt=1, transform=ccrs.PlateCarree(),norm=LogNorm()) 
cb = plt.colorbar(hb, ax=ax, orientation='vertical', pad=0.02, shrink=0.8)
cb.set_label('number of datapoints', fontsize=14)
gl=ax.gridlines(linewidth=0.2,color='grey',alpha=0.7,linestyle='-', draw_labels=True, x_inline= False,y_inline=False)
gl.xformatter=LONGITUDE_FORMATTER
gl.yformatter=LATITUDE_FORMATTER
gl.top_labels = False    # Disable top labels
gl.right_labels = False  # Disable right labels
ax.set_extent([-180, -40, 10, 80], crs=ccrs.PlateCarree())

ax.set_title('Spatial Data Density', fontsize=18, fontweight='bold')

plt.show()
_images/aad0e210b29cc76a3e1be9c92ba8ef192c72e781f24fd2c3fa86daa7c7abdafd.png