Quebrada Sonadora Nutrient Exports

Fluxes were calculated using the “loadflex” program (Appling et al. 2015), which is an enhancement of the widely used USGS “LOADEST” model (Runkel et al. 2004) as implemented in R. Both loadflex and LOADEST rely primarily on relationships between concentration of a solute and instantaneous discharge at the time of sampling to estimate concentrations when measured values of river chemistry are not available. In practice, and in our case, this often means that the 15-minute record of discharge at a given station is paired with an estimate of concentration based on weekly grab samples that span a range of discharge conditions.

First install some packages we will use.

In [ ]:
# install and load packages, you only need to do this once. 

rprofile_path = file.path(Sys.getenv("HOME"), ".Rprofile")
write('\noptions(repos=c(getOption(\'repos\'),
          CRAN=\'https://cloud.r-project.org\',
          USGS=\'https://owi.usgs.gov/R\'))\n',
      rprofile_path, 
      append =  TRUE)
install.packages('smwrBase')
install.packages('rloadest')
install.packages('car')
install.packages('unitted')
install.packages('dplyr')
packageurl <- "https://github.com/USGS-R/loadflex/archive/v1.0.1.tar.gz"
install.packages(packageurl, repos=NULL, type="source")
install.packages('dplyr')

load some r packages we will use.

In [ ]:
options(warn=-1)

library('loadflex')
library('rloadest')

You will need to include the data directory from this dataset in the same location you are running this notebook.

Initially here we will load discharge and weekly grab sample data from csv files. This notebook details nutrient fluxes for the Quebrada Sonadora sampling site which has a watershed area of 261.58 hectares.

In [1]:
QS_chem <- readRDS("data/QS_chem_UNH.rds")
QSDischargeShort <- readRDS("data/QS_Discharge_USGS.rds")
QSWatershedArea <- 261.5888

Calculate flux totals and mean annual concentration for Na.

In [2]:
QS_chemNa <- QS_chem[complete.cases(QS_chem["Na_mg_L"]),]
QS_chemNa <- QS_chemNa[complete.cases(QS_chemNa["CFS"]),]


meta <- metadata(constituent="Na_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1")


lr <- loadReg2(loadReg(Na_mg_L ~ model(1), data=QS_chemNa,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))


lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemNa, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

# preds_lrNa <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

preds_lcNa <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcNa <- aggregateSolute(preds_lcNa,meta, format="flux total", se.preds=preds_lcNa$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemNa, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lcNaConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcNaConc <- aggregateSolute(preds_lcNaConc,meta, format="conc", se.preds=preds_lcNaConc$se.pred, agg.by="calendar year")

aggs_lcNa$Na_Kg_Ha_yr <- aggs_lcNa$Flux_Total / QSWatershedArea
aggs_lcNa$Na_mg_l <- aggs_lcNaConc$Conc
FluxTotalsdf <- aggs_lcNa[, c('Na_mg_l','Na_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdf$Calendar_Year
FluxTotalsdf$Calendar_Year <- NULL
FluxTotalsdf <- round(FluxTotalsdf, 2)
FluxTotalsdf$Calendar_Year <- years

Calculate flux totals and mean annual concentration for Ca.

In [3]:
QS_chemCa <- QS_chem[complete.cases(QS_chem["Ca_mg_L"]),]
QS_chemCa <- QS_chemCa[complete.cases(QS_chemCa["CFS"]),]

meta <- metadata(constituent="Ca_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="Rio Piedras, PR")


lr <- loadReg2(loadReg(Ca_mg_L ~ model(1), data=QS_chemCa,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemCa, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lrCa <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrCa <- aggregateSolute(preds_lrCa,meta, format="flux total", se.preds=preds_lrCa$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemCa, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lcCaConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcCaConc <- aggregateSolute(preds_lcCaConc,meta, format="conc", se.preds=preds_lcCaConc$se.pred, agg.by="calendar year")

aggs_lrCa$Ca_Kg_Ha_yr <- aggs_lrCa$Flux_Total / QSWatershedArea
aggs_lrCa$Ca_mg_l <- aggs_lcCaConc$Conc
FluxTotalsdfCa <- aggs_lrCa[, c('Ca_mg_l','Ca_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfCa$Calendar_Year
FluxTotalsdfCa$Calendar_Year <- NULL
FluxTotalsdfCa <- round(FluxTotalsdfCa, 2)
FluxTotalsdfCa$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfCa,FluxTotalsdf, by = "Calendar_Year")
Warning message in loadReg(Ca_mg_L ~ model(1), data = QS_chemCa, flow = "CFS", dates = "date", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
In [4]:
# FluxTtotalsdf <- merge(FluxTtotalsdfCa,FluxTtotalsdf, by.x='Calendar_Year', by.y='Calendar_Year')
FluxTotals
Calendar_YearCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 2.70 56.914.76 109.37
2010 2.61 108.814.50 193.90
2011 2.83 129.894.66 225.30
2012 2.67 103.364.64 179.82
2013 2.52 93.104.60 171.38
2014 2.92 75.594.91 134.12
In [5]:
#Mg
QS_chemMg <- QS_chem[complete.cases(QS_chem["Mg_mg_L"]),]
QS_chemMg <- QS_chemMg[complete.cases(QS_chemMg["CFS"]),]

meta <- metadata(constituent="Mg_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(Mg_mg_L ~ model(1), data=QS_chemMg,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemMg, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lrMg <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrMg <- aggregateSolute(preds_lrMg,meta, format="flux total", se.preds=preds_lrMg$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemMg, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lcMgConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcMgConc <- aggregateSolute(preds_lcMgConc,meta, format="conc", se.preds=preds_lcMgConc$se.pred, agg.by="calendar year")

aggs_lrMg$Mg_Kg_Ha_yr <- aggs_lrMg$Flux_Total / QSWatershedArea
aggs_lrMg$Mg_mg_l <- aggs_lcMgConc$Conc
FluxTotalsdfMg <- aggs_lrMg[, c('Mg_mg_l','Mg_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfMg$Calendar_Year
FluxTotalsdfMg$Calendar_Year <- NULL
FluxTotalsdfMg <- round(FluxTotalsdfMg, 2)
FluxTotalsdfMg$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfMg,FluxTotals, by = "Calendar_Year")
FluxTotals
Warning message in loadReg(Mg_mg_L ~ model(1), data = QS_chemMg, flow = "CFS", dates = "date", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 1.69 34.81 2.70 56.914.76 109.37
2010 1.39 59.00 2.61 108.814.50 193.90
2011 1.40 62.26 2.83 129.894.66 225.30
2012 1.46 54.64 2.67 103.364.64 179.82
2013 1.58 54.93 2.52 93.104.60 171.38
2014 2.39 38.81 2.92 75.594.91 134.12
In [8]:
QS_chemNO3 <- QS_chem[complete.cases(QS_chem["NO3.N_ug_N_L"]),]
QS_chemNO3 <- QS_chemNO3[complete.cases(QS_chemNO3["CFS"]),]

meta <- metadata(constituent="NO3.N_ug_N_L", flow="CFS", 
                 dates="date", conc.units="ug/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(NO3.N_ug_N_L ~ model(1), data=QS_chemNO3,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="ug/L", load.units="kg"))
lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemNO3, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcNO3Conc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcNO3Conc <- aggregateSolute(preds_lcNO3Conc,meta, format="conc", se.preds=preds_lcNO3Conc$se.pred, agg.by="calendar year")

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemNO3, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrNO3 <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrNO3 <- aggregateSolute(preds_lrNO3,meta, format="flux total", se.preds=preds_lrNO3$se.pred, agg.by="calendar year")
aggs_lrNO3$NO3_Kg_Ha_yr <- aggs_lrNO3$Flux_Total / QSWatershedArea
aggs_lrNO3$NO3_mg_l <- aggs_lcNO3Conc$Conc / 1000.0
FluxTotalsdfNO3 <- aggs_lrNO3[, c('NO3_mg_l','NO3_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfNO3$Calendar_Year
FluxTotalsdfNO3$Calendar_Year <- NULL
FluxTotalsdfNO3 <- round(FluxTotalsdfNO3, 3)
FluxTotalsdfNO3$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfNO3,FluxTotals, by = "Calendar_Year")
FluxTotals
Warning message in loadReg(NO3.N_ug_N_L ~ model(1), data = QS_chemNO3, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearNO3_mg_l.xNO3_Kg_Ha_yr.xSO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 0.067 1.033 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 0.061 1.798 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 0.056 1.873 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 0.042 1.193 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 0.037 0.903 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.055 1.115 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [7]:
#SO4
QS_chemSO4 <- QS_chem[complete.cases(QS_chem["SO4_S_mg_L"]),]
QS_chemSO4 <- QS_chemSO4[complete.cases(QS_chemSO4["CFS"]),]

meta <- metadata(constituent="SO4_S_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(SO4_S_mg_L ~ model(1), data=QS_chemSO4,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))


lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemSO4, abs.or.rel.resids="absolute", interp.function=triangularInterpolation)

preds_lcSO4Conc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcSO4Conc <- aggregateSolute(preds_lcSO4Conc,meta, format="conc", se.preds=preds_lcSO4Conc$se.pred, agg.by="calendar year")

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemSO4, abs.or.rel.resids="relative", interp.function=triangularInterpolation)

preds_lrSO4 <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrSO4 <- aggregateSolute(preds_lrSO4,meta, format="flux total", se.preds=preds_lrSO4$se.pred, agg.by="calendar year")


aggs_lrSO4$SO4_Kg_Ha_yr <- aggs_lrSO4$Flux_Total / QSWatershedArea
aggs_lrSO4$SO4_mg_l <- aggs_lcSO4Conc$Conc
FluxTotalsdfSO4 <- aggs_lrSO4[, c('SO4_mg_l','SO4_Kg_Ha_yr','Calendar_Year')]


#round values 
years <- FluxTotalsdfSO4$Calendar_Year
FluxTotalsdfSO4$Calendar_Year <- NULL
FluxTotalsdfSO4 <- round(FluxTotalsdfSO4, 3)
FluxTotalsdfSO4$Calendar_Year <- years


FluxTotals <- merge(FluxTotalsdfSO4,FluxTotals, by = "Calendar_Year")
FluxTotals
Warning message in loadReg(SO4_S_mg_L ~ model(1), data = QS_chemSO4, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearSO4_mg_lSO4_Kg_Ha_yrNO3_mg_lNO3_Kg_Ha_yrMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [12]:
#F
QS_chemF <- QS_chem[complete.cases(QS_chem["Fluoride_mg_L"]),]
QS_chemF <- QS_chemF[complete.cases(QS_chemF["CFS"]),]
QSDischargeShort2<- QSDischargeShort[QSDischargeShort$date >= as.POSIXct("2013-01-01") 
                                                     & QSDischargeShort$date <= as.POSIXct("2014-12-31"),]

meta <- metadata(constituent="Fluoride_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(Fluoride_mg_L ~ model(1), data=QS_chemF,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))


lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemF, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcFConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort2, se.pred=TRUE, date=TRUE)

aggs_lcFConc <- aggregateSolute(preds_lcFConc,meta, format="conc", se.preds=preds_lcFConc$se.pred, agg.by="calendar year")

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemF, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrF <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort2, se.pred=TRUE, date=TRUE)

#get some really high wierd values just remove them. 
#preds_lrF <-preds_lrF[preds_lrF$fit < 50,]

aggs_lrF <- aggregateSolute(preds_lrF,meta, format="flux total", se.preds=preds_lrF$se.pred, agg.by="calendar year")
aggs_lcFConc <- aggregateSolute(preds_lcFConc,meta, format="conc", se.preds=preds_lcFConc$se.pred, agg.by="calendar year")
aggs_lrF$F_Kg_Ha_yr <- aggs_lrF$Flux_Total / QSWatershedArea
aggs_lrF$F_mg_l <- aggs_lcFConc$Conc
FluxTotalsdfF <- aggs_lrF[, c('F_mg_l','F_Kg_Ha_yr','Calendar_Year')]

# round values
years <- FluxTotalsdfF$Calendar_Year
FluxTotalsdfF$Calendar_Year <- NULL
FluxTotalsdfF <- round(FluxTotalsdfF, 3)
FluxTotalsdfF$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfF,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Calendar_YearF_mg_lF_Kg_Ha_yrNO3_mg_lNO3_Kg_Ha_yrSO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2013 0.017 0.589 0.037 0.903 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.019 0.452 0.055 1.115 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
2009 NA NA 0.067 1.033 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 NA NA 0.061 1.798 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 NA NA 0.056 1.873 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 NA NA 0.042 1.193 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
In [ ]:
QS_chemSiO2 <- QS_chem[complete.cases(QS_chem["SiO2_mg_L"]),]
QS_chemSiO2 <- QS_chemSiO2[complete.cases(QS_chemSiO2["CFS"]),]

meta <- metadata(constituent="SiO2_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(SiO2_mg_L ~ model(1), data=QS_chemSiO2,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemSiO2, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lcSiO2Conc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcSiO2Conc <- aggregateSolute(preds_lcSiO2Conc,meta, format="conc", se.preds=preds_lcSiO2Conc$se.pred, agg.by="calendar year")




lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemSiO2, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lrSiO2 <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrSiO2 <- aggregateSolute(preds_lrSiO2,meta, format="flux total", se.preds=preds_lrSiO2$se.pred, agg.by="calendar year")

aggs_lrSiO2$SiO2_Kg_Ha_yr <- aggs_lrSiO2$Flux_Total / QSWatershedArea
aggs_lrSiO2$SiO2_mg_l <- aggs_lcSiO2Conc$Conc
FluxTotalsdfSiO2 <- aggs_lrSiO2[, c('SiO2_mg_l','SiO2_Kg_Ha_yr','Calendar_Year')]


#round values 
years <- FluxTotalsdfSiO2$Calendar_Year
FluxTotalsdfSiO2$Calendar_Year <- NULL
FluxTotalsdfSiO2 <- round(FluxTotalsdfSiO2, 2)
FluxTotalsdfSiO2$Calendar_Year <- years


FluxTotals <- merge(FluxTotalsdfSiO2,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
In [13]:
#PO4
QS_chemPO4 <- QS_chem[complete.cases(QS_chem["PO4.P_ug_P_L"]),]
QS_chemPO4 <- QS_chemPO4[complete.cases(QS_chemPO4["CFS"]),]

meta <- metadata(constituent="PO4.P_ug_P_L", flow="CFS", 
                 dates="date", conc.units="ug/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(PO4.P_ug_P_L ~ model(1), data=QS_chemPO4,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="ug/L", load.units="kg"))
lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemPO4, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcPO4Conc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcPO4Conc <- aggregateSolute(preds_lcPO4Conc,meta, format="conc", se.preds=preds_lcPO4Conc$se.pred, agg.by="calendar year")

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemPO4, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrPO4 <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrPO4 <- aggregateSolute(preds_lrPO4,meta, format="flux total", se.preds=preds_lrPO4$se.pred, agg.by="calendar year")

aggs_lrPO4$PO4_Kg_Ha_yr <- aggs_lrPO4$Flux_Total / QSWatershedArea
aggs_lrPO4$PO4_mg_l <- aggs_lcPO4Conc$Conc / 1000.0
FluxTotalsdfPO4 <- aggs_lrPO4[, c('PO4_mg_l','PO4_Kg_Ha_yr','Calendar_Year')]


#round values 
years <- FluxTotalsdfPO4$Calendar_Year
FluxTotalsdfPO4$Calendar_Year <- NULL
FluxTotalsdfPO4$PO4_mg_l <- round(FluxTotalsdfPO4$PO4_mg_l, 3)
FluxTotalsdfPO4$PO4_Kg_Ha_yr <- round(FluxTotalsdfPO4$PO4_Kg_Ha_yr, 2)

FluxTotalsdfPO4$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfPO4,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(PO4.P_ug_P_L ~ model(1), data = QS_chemPO4, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearPO4_mg_lPO4_Kg_Ha_yrF_mg_lF_Kg_Ha_yrNO3_mg_lNO3_Kg_Ha_yrSO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 0.003 0.06 NA NA 0.067 1.033 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 0.002 0.11 NA NA 0.061 1.798 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 0.002 0.12 NA NA 0.056 1.873 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 0.002 0.09 NA NA 0.042 1.193 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 0.002 0.10 0.017 0.589 0.037 0.903 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.003 0.09 0.019 0.452 0.055 1.115 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [14]:
#NH4 

QS_chemNH4 <- QS_chem[complete.cases(QS_chem["NH4.N_ug_L"]),]
QS_chemNH4 <- QS_chemNH4[complete.cases(QS_chemNH4["CFS"]),]

meta <- metadata(constituent="NH4.N_ug_L", flow="CFS", 
                 dates="date", conc.units="ug/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(NH4.N_ug_L ~ model(1), data=QS_chemNH4,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="ug/L", load.units="kg"))
lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemNH4, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcNH4Conc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcNH4Conc <- aggregateSolute(preds_lcNH4Conc,meta, format="conc", se.preds=preds_lcNH4Conc$se.pred, agg.by="calendar year")

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemNH4, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrNH4 <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrNH4 <- aggregateSolute(preds_lrNH4,meta, format="flux total", se.preds=preds_lrNH4$se.pred, agg.by="calendar year")

aggs_lrNH4$NH4_Kg_Ha_yr <- aggs_lrNH4$Flux_Total / QSWatershedArea
aggs_lrNH4$NH4_mg_l <- aggs_lcNH4Conc$Conc
FluxTotalsdfNH4 <- aggs_lrNH4[, c('NH4_mg_l','NH4_Kg_Ha_yr','Calendar_Year')]


#round values 
years <- FluxTotalsdfNH4$Calendar_Year
FluxTotalsdfNH4$Calendar_Year <- NULL
FluxTotalsdfNH4$NH4_mg_l <- round(FluxTotalsdfNH4$NH4_mg_l, 3)
FluxTotalsdfNH4$NH4_Kg_Ha_yr <- round(FluxTotalsdfNH4$NH4_Kg_Ha_yr, 2)
FluxTotalsdfNH4$Calendar_Year <- years


FluxTotals <- merge(FluxTotalsdfNH4,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(NH4.N_ug_L ~ model(1), data = QS_chemNH4, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearNH4_mg_lNH4_Kg_Ha_yrPO4_mg_lPO4_Kg_Ha_yrF_mg_lF_Kg_Ha_yrNO3_mg_lNO3_Kg_Ha_yrSO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 5.801 0.17 0.003 0.06 NA NA 0.067 1.033 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 8.741 0.43 0.002 0.11 NA NA 0.061 1.798 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 4.176 0.22 0.002 0.12 NA NA 0.056 1.873 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 4.791 0.17 0.002 0.09 NA NA 0.042 1.193 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 5.029 0.19 0.002 0.10 0.017 0.589 0.037 0.903 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 5.423 0.19 0.003 0.09 0.019 0.452 0.055 1.115 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [15]:
QS_chemTDN <- QS_chem[complete.cases(QS_chem["TDN_mg_L"]),]
QS_chemTDN <- QS_chemTDN[complete.cases(QS_chemTDN["CFS"]),]

meta <- metadata(constituent="TDN_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(TDN_mg_L ~ model(1), data=QS_chemTDN,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemTDN, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrTDN <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrTDN <- aggregateSolute(preds_lrTDN,meta, format="flux total", se.preds=preds_lrTDN$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemTDN, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcTDNConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcTDNConc <- aggregateSolute(preds_lcTDNConc,meta, format="conc", se.preds=preds_lcTDNConc$se.pred, agg.by="calendar year")

aggs_lrTDN$TDN_Kg_Ha_yr <- aggs_lrTDN$Flux_Total / QSWatershedArea
aggs_lrTDN$TDN_mg_l <- aggs_lcTDNConc$Conc
FluxTotalsdfTDN <- aggs_lrTDN[, c('TDN_mg_l','TDN_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfTDN$Calendar_Year
FluxTotalsdfTDN$Calendar_Year <- NULL
FluxTotalsdfTDN <- round(FluxTotalsdfTDN, 2)
FluxTotalsdfTDN$Calendar_Year <- years


FluxTotals <- merge(FluxTotalsdfTDN,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(TDN_mg_L ~ model(1), data = QS_chemTDN, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearTDN_mg_lTDN_Kg_Ha_yrNH4_mg_lNH4_Kg_Ha_yrPO4_mg_lPO4_Kg_Ha_yrF_mg_lF_Kg_Ha_yrNO3_mg_l...SO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 0.14 3.16 5.801 0.17 0.003 0.06 NA NA 0.067 ... 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 0.18 6.38 8.741 0.43 0.002 0.11 NA NA 0.061 ... 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 0.14 6.53 4.176 0.22 0.002 0.12 NA NA 0.056 ... 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 0.09 3.27 4.791 0.17 0.002 0.09 NA NA 0.042 ... 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 0.10 3.47 5.029 0.19 0.002 0.10 0.017 0.589 0.037 ... 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.12 3.36 5.423 0.19 0.003 0.09 0.019 0.452 0.055 ... 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [16]:
#DON
QS_chemDON <- QS_chem[complete.cases(QS_chem["DON_mg_L"]),]
QS_chemDON <- QS_chemDON[complete.cases(QS_chemDON["CFS"]),]

meta <- metadata(constituent="DON_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(DON_mg_L ~ model(1), data=QS_chemDON,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemDON, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrDON <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrDON <- aggregateSolute(preds_lrDON,meta, format="flux total", se.preds=preds_lrDON$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemDON, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcDONConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcDONConc <- aggregateSolute(preds_lcDONConc,meta, format="conc", se.preds=preds_lcDONConc$se.pred, agg.by="calendar year")


aggs_lrDON$DON_Kg_Ha_yr <- aggs_lrDON$Flux_Total / QSWatershedArea
aggs_lrDON$DON_mg_l <- aggs_lcDONConc$Conc
FluxTotalsdfDON <- aggs_lrDON[, c('DON_mg_l','DON_Kg_Ha_yr','Calendar_Year')]


#round values 
years <- FluxTotalsdfDON$Calendar_Year
FluxTotalsdfDON$Calendar_Year <- NULL
FluxTotalsdfDON <- round(FluxTotalsdfDON, 2)
FluxTotalsdfDON$Calendar_Year <- years


FluxTotals <- merge(FluxTotalsdfDON,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(DON_mg_L ~ model(1), data = QS_chemDON, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearDON_mg_lDON_Kg_Ha_yrTDN_mg_lTDN_Kg_Ha_yrNH4_mg_lNH4_Kg_Ha_yrPO4_mg_lPO4_Kg_Ha_yrF_mg_l...SO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 0.07 1.99 0.14 3.16 5.801 0.17 0.003 0.06 NA ... 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 0.06 2.93 0.18 6.38 8.741 0.43 0.002 0.11 NA ... 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 0.07 4.47 0.14 6.53 4.176 0.22 0.002 0.12 NA ... 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 0.04 1.91 0.09 3.27 4.791 0.17 0.002 0.09 NA ... 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 0.06 2.50 0.10 3.47 5.029 0.19 0.002 0.10 0.017 ... 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.06 1.91 0.12 3.36 5.423 0.19 0.003 0.09 0.019 ... 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [17]:
#DOC
QS_chemDOC <- QS_chem[complete.cases(QS_chem["DOC_mg_L"]),]
QS_chemDOC <- QS_chemDOC[complete.cases(QS_chemDOC["CFS"]),]

meta <- metadata(constituent="DOC_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(DOC_mg_L ~ model(1), data=QS_chemDOC,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemDOC, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lrDOC <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrDOC <- aggregateSolute(preds_lrDOC,meta, format="flux total", se.preds=preds_lrDOC$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemDOC, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lcDOCConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcDOCConc <- aggregateSolute(preds_lcDOCConc,meta, format="conc", se.preds=preds_lcDOCConc$se.pred, agg.by="calendar year")

aggs_lrDOC$DOC_Kg_Ha_yr <- aggs_lrDOC$Flux_Total / QSWatershedArea
aggs_lrDOC$DOC_mg_l <- aggs_lcDOCConc$Conc
FluxTotalsdfDOC <- aggs_lrDOC[, c('DOC_mg_l','DOC_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfDOC$Calendar_Year
FluxTotalsdfDOC$Calendar_Year <- NULL
FluxTotalsdfDOC <- round(FluxTotalsdfDOC, 2)
FluxTotalsdfDOC$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfDOC,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(DOC_mg_L ~ model(1), data = QS_chemDOC, flow = "CFS", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearDOC_mg_lDOC_Kg_Ha_yrDON_mg_lDON_Kg_Ha_yrTDN_mg_lTDN_Kg_Ha_yrNH4_mg_lNH4_Kg_Ha_yrPO4_mg_l...SO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 2.22 98.320.07 1.99 0.14 3.16 5.801 0.17 0.003 ... 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 2.05 159.900.06 2.93 0.18 6.38 8.741 0.43 0.002 ... 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 1.95 152.430.07 4.47 0.14 6.53 4.176 0.22 0.002 ... 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 1.74 97.330.04 1.91 0.09 3.27 4.791 0.17 0.002 ... 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 1.76 143.950.06 2.50 0.10 3.47 5.029 0.19 0.002 ... 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 1.61 62.290.06 1.91 0.12 3.36 5.423 0.19 0.003 ... 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [18]:
#Cl
QS_chemCl <- QS_chem[complete.cases(QS_chem["Cl_mg_L"]),]
QS_chemCl <- QS_chemCl[complete.cases(QS_chemCl["CFS"]),]

meta <- metadata(constituent="Cl_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(Cl_mg_L ~ model(1), data=QS_chemCl,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemCl, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lrCl <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrCl <- aggregateSolute(preds_lrCl,meta, format="flux total", se.preds=preds_lrCl$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemCl, abs.or.rel.resids="relative", interp.function=rectangularInterpolation)

preds_lcClConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lcClConc <- aggregateSolute(preds_lcClConc,meta, format="conc", se.preds=preds_lcClConc$se.pred, agg.by="calendar year")


aggs_lrCl$Cl_Kg_Ha_yr <- aggs_lrCl$Flux_Total / QSWatershedArea
aggs_lrCl$Cl_mg_l <- aggs_lcClConc$Conc
FluxTotalsdfCl <- aggs_lrCl[, c('Cl_mg_l','Cl_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfCl$Calendar_Year
FluxTotalsdfCl$Calendar_Year <- NULL
FluxTotalsdfCl <- round(FluxTotalsdfCl, 2)
FluxTotalsdfCl$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfCl,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(Cl_mg_L ~ model(1), data = QS_chemCl, flow = "CFS", dates = "date", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearCl_mg_lCl_Kg_Ha_yrDOC_mg_lDOC_Kg_Ha_yrDON_mg_lDON_Kg_Ha_yrTDN_mg_lTDN_Kg_Ha_yrNH4_mg_l...SO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 7.22 181.602.22 98.320.07 1.99 0.14 3.16 5.801 ... 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 6.69 306.032.05 159.900.06 2.93 0.18 6.38 8.741 ... 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 6.72 327.981.95 152.430.07 4.47 0.14 6.53 4.176 ... 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 7.02 281.041.74 97.330.04 1.91 0.09 3.27 4.791 ... 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 7.25 290.831.76 143.950.06 2.50 0.10 3.47 5.029 ... 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 7.38 213.481.61 62.290.06 1.91 0.12 3.36 5.423 ... 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12
In [19]:
QS_chemK <- QS_chem[complete.cases(QS_chem["K_mg_L"]),]
QS_chemK <- QS_chemK[complete.cases(QS_chemK["CFS"]),]

meta <- metadata(constituent="K_mg_L", flow="CFS", 
                 dates="date", conc.units="mg/L", flow.units="cfs", load.units="kg", 
                 load.rate.units="kg d^-1", station="QS, PR")


lr <- loadReg2(loadReg(K_mg_L ~ model(1), data=QS_chemK,
                       flow="CFS", dates="date", time.step="instantaneous", 
                       flow.units="cfs", conc.units="mg/L", load.units="kg"))

lc <- loadComp(reg.model=lr, interp.format="flux", 
               interp.data=QS_chemK, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lrK <- predictSolute(lc, "flux",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)

aggs_lrK <- aggregateSolute(preds_lrK,meta, format="flux total", se.preds=preds_lrK$se.pred, agg.by="calendar year")

lcConc <- loadComp(reg.model=lr, interp.format="conc", 
                   interp.data=QS_chemK, abs.or.rel.resids="absolute", interp.function=rectangularInterpolation)

preds_lcKConc <- predictSolute(lcConc, "conc",interval="prediction", QSDischargeShort, se.pred=TRUE, date=TRUE)


aggs_lrK$K_Kg_Ha_yr <- aggs_lrK$Flux_Total / QSWatershedArea
aggs_lcKConc <- aggregateSolute(preds_lcKConc,meta, format="conc", se.preds=preds_lcKConc$se.pred, agg.by="calendar year")
aggs_lrK$K_mg_l <- aggs_lcKConc$Conc
FluxTotalsdfK <- aggs_lrK[, c('K_mg_l','K_Kg_Ha_yr','Calendar_Year')]

#round values 
years <- FluxTotalsdfK$Calendar_Year
FluxTotalsdfK$Calendar_Year <- NULL
FluxTotalsdfK <- round(FluxTotalsdfK, 2)
FluxTotalsdfK$Calendar_Year <- years

FluxTotals <- merge(FluxTotalsdfK,FluxTotals, by = "Calendar_Year", all.y=TRUE)
FluxTotals
Warning message in loadReg(K_mg_L ~ model(1), data = QS_chemK, flow = "CFS", dates = "date", :
"The minimum spacing between observed loads is 6 days. The time between observations should be at least  7 days to avoid autocorrelation issues."
Calendar_YearK_mg_lK_Kg_Ha_yrCl_mg_lCl_Kg_Ha_yrDOC_mg_lDOC_Kg_Ha_yrDON_mg_lDON_Kg_Ha_yrTDN_mg_l...SO4_mg_lSO4_Kg_Ha_yrNO3_mg_l.yNO3_Kg_Ha_yr.yMg_mg_lMg_Kg_Ha_yrCa_mg_lCa_Kg_Ha_yrNa_mg_lNa_Kg_Ha_yr
2009 0.29 6.99 7.22 181.602.22 98.320.07 1.99 0.14 ... 0.55 15.21 0.07 1.03 1.69 34.81 2.70 56.914.76 109.37
2010 0.27 11.45 6.69 306.032.05 159.900.06 2.93 0.18 ... 0.53 25.03 0.06 1.80 1.39 59.00 2.61 108.814.50 193.90
2011 0.21 10.49 6.72 327.981.95 152.430.07 4.47 0.14 ... 0.52 25.99 0.06 1.87 1.40 62.26 2.83 129.894.66 225.30
2012 0.24 9.55 7.02 281.041.74 97.330.04 1.91 0.09 ... 0.52 20.92 0.04 1.19 1.46 54.64 2.67 103.364.64 179.82
2013 0.23 8.80 7.25 290.831.76 143.950.06 2.50 0.10 ... 0.51 20.67 0.04 0.90 1.58 54.93 2.52 93.104.60 171.38
2014 0.30 8.55 7.38 213.481.61 62.290.06 1.91 0.12 ... 0.54 15.46 0.05 1.11 2.39 38.81 2.92 75.594.91 134.12

Format the table of fluxes and concentrations.

In [20]:
rownames(FluxTotals) <- FluxTotals$Calendar_Year
FluxTotals2 <- as.data.frame(t(FluxTotals))
FluxTotals2$Calendar_Year <- NULL
FluxTotals2
200920102011201220132014
Calendar_Year2009 2010 2011 2012 2013 2014
K_mg_l0.29 0.27 0.21 0.24 0.23 0.30
K_Kg_Ha_yr 6.99 11.45 10.49 9.55 8.80 8.55
Cl_mg_l7.22 6.69 6.72 7.02 7.25 7.38
Cl_Kg_Ha_yr181.60306.03327.98281.04290.83213.48
DOC_mg_l2.22 2.05 1.95 1.74 1.76 1.61
DOC_Kg_Ha_yr 98.32159.90152.43 97.33143.95 62.29
DON_mg_l0.07 0.06 0.07 0.04 0.06 0.06
DON_Kg_Ha_yr1.99 2.93 4.47 1.91 2.50 1.91
TDN_mg_l0.14 0.18 0.14 0.09 0.10 0.12
TDN_Kg_Ha_yr3.16 6.38 6.53 3.27 3.47 3.36
NH4_mg_l5.801 8.741 4.176 4.791 5.029 5.423
NH4_Kg_Ha_yr0.17 0.43 0.22 0.17 0.19 0.19
PO4_mg_l0.003 0.002 0.002 0.002 0.002 0.003
PO4_Kg_Ha_yr0.06 0.11 0.12 0.09 0.10 0.09
F_mg_lNA NA NA NA 0.017 0.019
F_Kg_Ha_yrNA NA NA NA 0.589 0.452
NO3_mg_l0.067 0.061 0.056 0.042 0.037 0.055
NO3_Kg_Ha_yr1.033 1.798 1.873 1.193 0.903 1.115
SO4_mg_l0.55 0.53 0.52 0.52 0.51 0.54
SO4_Kg_Ha_yr15.21 25.03 25.99 20.92 20.67 15.46
NO3_mg_l.y0.07 0.06 0.06 0.04 0.04 0.05
NO3_Kg_Ha_yr.y1.03 1.80 1.87 1.19 0.90 1.11
Mg_mg_l1.69 1.39 1.40 1.46 1.58 2.39
Mg_Kg_Ha_yr34.81 59.00 62.26 54.64 54.93 38.81
Ca_mg_l2.70 2.61 2.83 2.67 2.52 2.92
Ca_Kg_Ha_yr 56.91108.81129.89103.36 93.10 75.59
Na_mg_l4.76 4.50 4.66 4.64 4.60 4.91
Na_Kg_Ha_yr109.37193.90225.30179.82171.38134.12
In [ ]: