use nixmass R package from python

calculate_deltasnow.R

# this file can be sourced from python with rpy2

library(nixmass)

calculate_deltasnow <- function(df_in){
  #' Calculate nixmass delta snow and return dataframe instead of nixmass object
  #'
  #' @param df_in The dataframe containing the data 
  #
  if(df_in$hs[1]!=0) {  
    # if first HS datapoint is not zero, raise second datapoint by the half of first and set first to zero.
    df_in$hs[2]<-df_in$hs[2]+(df_in$hs[1]/2)
    df_in$hs[1]<-0}
  
  # call nixmass function with delta.snow option
  SWEnixmass<-nixmass(df_in, model="delta.snow", verbose=F)
  
  # create a data.frame containing only SWE from the nixmass object.
  SWE_delta_snow<-SWEnixmass$swe$delta.snow
  date <- SWEnixmass$date
  df_out <- data.frame(date, SWE_delta_snow)
  return(df_out)
}

nixmass_from_python.py

# -*- coding: utf-8 -*-
"""
Call nixmass R package from python with rpy2.

Sample HS data has been downloaded from the nixmass cran mirror:
https://github.com/cran/nixmass/blob/master/data/hsdata.RData

Feel free copy if something is helpful for your usecase.

Works for me with the following versions:
python: '3.7.6'
R: '3.6.1'
rpy2: '3.4.5'
pandas: '1.1.3'

@author: Johannes Aschauer
"""
import pandas as pd

import rpy2.robjects as ro
import rpy2.robjects.packages as rpackages
from rpy2.robjects.vectors import StrVector
from rpy2.robjects import pandas2ri

# R package setup and installaton:
REQUIRED_R_PACKAGES = ['nixmass']
utils = rpackages.importr('utils')
utils.chooseCRANmirror(ind=1)
PACKAGES_TO_INSTALL = [p for p in REQUIRED_R_PACKAGES if not rpackages.isinstalled(p)]
if len(PACKAGES_TO_INSTALL) > 0:
    utils.install_packages(StrVector(PACKAGES_TO_INSTALL))

# activate automatic pandas to R conversion
pandas2ri.activate()

# Source the modified deltasnow function from R
ro.r['source']('calculate_deltasnow.R')
# create a pointer to a python function
calculate_deltasnow = ro.globalenv['calculate_deltasnow']


def append_delta_snow_SWE_column(df, hs_unit='m', swe_unit='mm'):
    """
    Calculate SWE from HS with nixmass R package and append it to input.

    Parameters
    ----------
    df : pd.DataFrame
        input data needs to have a column 'HS' or 'hs' and pd.DatetimeIndex.

    hs_unit : str of {'mm', 'cm', 'm'}
        Input HS unit. Will be transformed to [m] for the nixmass call. The
        default is 'm'

    swe_unit : str of {'mm', 'cm', 'm'}
        Output unit of the appended SWE column

    Returns
    -------
    df_modified : pd.DataFrame
        includes column 'SWE_delta_snow' in the defined output `swe_unit`

    """
    df = df.copy()
    assert hs_unit in {'mm', 'cm', 'm'}
    assert swe_unit in {'mm', 'cm', 'm'}
    assert 'hs' in [c.lower() for c in df.columns]
    assert isinstance(df.index, pd.DatetimeIndex)
    
    unit_factor = {'mm': 0.001,
                   'cm': 0.01,
                   'm': 1.0}
    
    # create water year column for grouping
    df['wy'] = df.index.year.where(df.index.month<9, df.index.year+1)
    to_concatenate = []
    for wy, data in df.groupby('wy'):
        # nixmass needs columns 'date' and 'hs'
        df_R_in = (data
                 .reset_index(drop=False)
                 .rename(columns={'index': 'date',
                                  data.index.name: 'date',
                                  'HS': 'hs'},
                         errors='ignore')
                 .loc[:,['date', 'hs']])
        # input HS for nixmass needs to be in [m]
        df_R_in['hs'] = df_R_in['hs']*unit_factor[hs_unit]
        # nixmass needs date as strings
        df_R_in['date'] = df_R_in.date.astype(str)
        # call R function
        df_R_out = calculate_deltasnow(df_R_in)
        # reconvert datetime format
        df_R_out['date'] = pd.to_datetime(df_R_out['date'])
        df_R_out = df_R_out.set_index('date')
        # append modeled SWE to data
        data['SWE_delta_snow'] = df_R_out['SWE_delta_snow']
        # nixmass returns SWE in [mm]
        data['SWE_delta_snow'] = data['SWE_delta_snow']*0.001/unit_factor[swe_unit]
        to_concatenate.append(data)

    df_modified = (pd.concat(to_concatenate, axis=0)
                   .drop(columns=['wy'])
                   )
    return df_modified


if __name__ == '__main__':
    # Test demo:
    # load sample data for testing (HS given in [m])
    data = pd.read_csv("sample_hsdata.csv",
                       index_col='date',
                       parse_dates=['date'])

    data['hs'] = data['hs']*100 # many times HS comes in cm.
    data_out = append_delta_snow_SWE_column(data, hs_unit='cm', swe_unit='mm')
    data_out.plot()
    
    # the `append_delta_snow_SWE_column` function can also very conveniently be used in .pipe() 
    # calls e.g. before monthly mean resampling:
    data2 = (pd.read_csv("sample_hsdata.csv",
                        index_col='date',
                        parse_dates=['date'])
            .pipe(append_delta_snow_SWE_column, hs_unit='m', swe_unit='m')
            .resample('M').mean()
            )
    data2.plot()
    
    

sample_hsdata.csv