doi-usgs / egretci Goto Github PK
View Code? Open in Web Editor NEWA bootstrap method for estimating uncertainty of water quality trends. https://doi-usgs.github.io/EGRETci/
Home Page: http://doi-usgs.github.io/EGRETci/
License: Other
A bootstrap method for estimating uncertainty of water quality trends. https://doi-usgs.github.io/EGRETci/
Home Page: http://doi-usgs.github.io/EGRETci/
License: Other
The nBoot runs are independent, I think it would be OK. The screen messages would be either impossible or harder to produce, but that's probably a small price to pay. Adding the complexity is the real question.
Hello,
Getting this reproducible error while running wBT and runPairs:
Error in if (z) "Reject Ho" else "Do Not Reject Ho" : missing value where TRUE/FALSE needed
Attached are example input text files used in model runs that produce error and my code as R script.
When using either input text file "as is", I get the error and change calculations are significantly off (Flux Change is actually NA):
"Susquehanna River at Danville, PA Nitrate-Nitrite
Water Year
Bootstrap process, for change from Water Year1986 to Water Year2019
data set runs from 1986 to Year 2019
Bootstrap block length in days 200
bootBreak is 3 confStop is 0.7
WRTDS estimated concentration change is 25256 mg/L
WRTDS estimated flux change is NA 10^6 kg/yr"
When I inputted "fake" results that complete the 1986 water year (in my case, I actually added fake sample results just for September or October), I no longer get the error.
Reproduced this same error for other input files and other sites.
Must the input sample result dates exactly coincide with the period of analysis (PA) in order for wBT to run without this error?
Thank you!
01540500_trend_git.zip
01540500dn - Original.txt
01540500dn.txt
I am currently rerunning a model that looks at a 10 yr moving window with nBoot = 200, and I keep noticing odd iteration values. Here is some output from 1998-2007:
iBoot, xConc and xFlux 127 -8.615097e+137 -2.975343e+136
iBoot, xConc and xFlux 128 -0.05028564 -1.459841
iBoot, xConc and xFlux 129 0.03666643 -0.1706541
iBoot, xConc and xFlux 130 0.05564319 0.03039397
iBoot, xConc and xFlux 131 0.04101899 -0.1198773
iBoot, xConc and xFlux 132 0.04211883 0.1540202
iBoot, xConc and xFlux 133 -0.09984916 -2.680429
iBoot, xConc and xFlux 134 -0.007287396 -0.3895768
iBoot, xConc and xFlux 135 0.002760677 -0.7174396
iBoot, xConc and xFlux 136 -0.001037794 -0.5118047
iBoot, xConc and xFlux 137 0.0370652 -0.2164727
iBoot, xConc and xFlux 138 -0.02083882 -0.5030731
iBoot, xConc and xFlux 139 0.02780439 -0.294982
iBoot, xConc and xFlux 140 0.003192806 -0.4554299
iBoot, xConc and xFlux 141 -0.05443917 -1.415746
iBoot, xConc and xFlux 142 0.03259536 -0.466467
iBoot, xConc and xFlux 143 0.01166099 -0.1782373
iBoot, xConc and xFlux 144 0.03090554 -0.2704868
iBoot, xConc and xFlux 145 -0.01288284 -1.037976
iBoot, xConc and xFlux 146 0.1241648 0.8563349
iBoot, xConc and xFlux 147 -0.03576329 -1.28158
iBoot, xConc and xFlux 148 2.897706e+30 1.57178e+32
iBoot, xConc and xFlux 149 0.08204259 0.3758117
iBoot, xConc and xFlux 150 -0.01351132 -0.4975678
iBoot, xConc and xFlux 151 0.01373099 -0.1204261
iBoot, xConc and xFlux 152 0.0691038 0.1389613
iBoot, xConc and xFlux 153 0.009624788 -0.4152171
iBoot, xConc and xFlux 154 -0.06310256 -1.229834
iBoot, xConc and xFlux 155 -0.03640907 -1.312824
iBoot, xConc and xFlux 156 0.02100833 -0.2485315
iBoot, xConc and xFlux 157 0.01610422 -0.3146633
iBoot, xConc and xFlux 158 -0.01941998 -1.072092
iBoot, xConc and xFlux 159 -0.06849002 -1.519826
iBoot, xConc and xFlux 160 0.06798843 0.2361849
iBoot, xConc and xFlux 161 0.01927036 -0.3341333
iBoot, xConc and xFlux 162 0.005404621 -0.6666091
iBoot, xConc and xFlux 163 0.04532531 -0.4035346
iBoot, xConc and xFlux 164 0.1073366 0.2855145
iBoot, xConc and xFlux 165 -0.02089126 -1.265881
iBoot, xConc and xFlux 166 -0.03983874 -1.183583
iBoot, xConc and xFlux 167 -1.307639e+79 -4.516112e+77
iBoot, xConc and xFlux 168 0.02232687 -0.2135413
iBoot, xConc and xFlux 169 -0.02690708 -1.169837
iBoot, xConc and xFlux 170 0.05313271 0.02626181
iBoot, xConc and xFlux 171 0.06698775 0.1982021
iBoot, xConc and xFlux 172 0.03341398 -0.7828989
iBoot, xConc and xFlux 173 0.07338573 0.1200608
iBoot, xConc and xFlux 174 0.07082601 0.122985
iBoot, xConc and xFlux 175 -7.434473e+101 -2.567598e+100
iBoot, xConc and xFlux 176 0.07073133 0.1410189
iBoot, xConc and xFlux 177 0.01842965 -0.3556271
iBoot, xConc and xFlux 178 -0.09302633 -2.252569
iBoot, xConc and xFlux 179 0.07934037 0.3939216
iBoot, xConc and xFlux 180 0.01684882 -0.2692689
iBoot, xConc and xFlux 181 0.02964044 -0.3813411
iBoot, xConc and xFlux 182 0.004074809 -0.3463523
iBoot, xConc and xFlux 183 -0.05169548 -1.534628
iBoot, xConc and xFlux 184 0.03390354 -0.1392435
iBoot, xConc and xFlux 185 0.06523865 0.07551456
iBoot, xConc and xFlux 186 0.0444889 -0.09759831
iBoot, xConc and xFlux 187 0.006336869 -0.4673323
iBoot, xConc and xFlux 188 -0.02044203 -1.04437
iBoot, xConc and xFlux 189 0.03901774 0.08431938
iBoot, xConc and xFlux 190 -0.04834053 -1.364399
iBoot, xConc and xFlux 191 -0.05268685 -1.252619
iBoot, xConc and xFlux 192 0.07842773 0.221899
iBoot, xConc and xFlux 193 -0.02316459 -1.329317
iBoot, xConc and xFlux 194 -0.09582452 -1.921
iBoot, xConc and xFlux 195 0.1049984 0.6393139
iBoot, xConc and xFlux 196 0.02325422 -0.2404297
iBoot, xConc and xFlux 197 -0.06165737 -1.684607
iBoot, xConc and xFlux 198 0.0514432 -0.08819077
iBoot, xConc and xFlux 199 -0.04280059 -1.591571
iBoot, xConc and xFlux 200 0.004167874 -0.5705411
It took 203 iterations to achieve 200 sucessful runs.
Maumee River at Waterville OH
TKN through 2007
Water Year
Change estimates are for 2007 minus 1998
Should we reject Ho that Flow Normalized Concentration Trend = 0 ? Do Not Reject Ho
best estimate of change in concentration is 0.0194 mg/L
Lower and Upper 90% CIs -1.19e-01 7.84e-02
also 95% CIs -1.46e+101 1.04e-01
and 50% CIs -2.75e-02 3.97e-02
approximate two-sided p-value for Conc 0.82
Likelihood that Flow Normalized Concentration is trending up = 0.59 is trending down = 0.41
Should we reject Ho that Flow Normalized Flux Trend = 0 ? Do Not Reject Ho
best estimate of change in flux is -0.409 10^6 kg/year
Lower and Upper 90% CIs -2.66e+00 2.80e-01
also 95% CIs -5.04e+99 4.97e-01
and 50% CIs -1.13e+00 -1.22e-01
approximate two-sided p-value for Flux 0.36
Likelihood that Flow Normalized Flux is trending up = 0.177 is trending down = 0.823
Do you have any idea why some of the iterations are so off? Thank you!
dif <- flexBoot$z11[1:500] - flexBoot$z00[1:500]
summary(dif)
sortDif <- sort(dif)
trueCI <- quantile(sortDif, probs=c(0.05,0.95), type = 6)
trueWidth <- trueCI[2] - trueCI[1]
nReps <- 500
bootWidth <- rep(NA,nReps)
for (iRep in 1:nReps){
difB <- sample(dif,500,replace = TRUE)
difB <- sort(difB)
CI <- quantile(difB, probs=c(0.05, 0.95), type = 6)
bootWidth[iRep] <- CI[2] - CI[1]
}
print(sd(bootWidth))
summary(bootWidth)
cat("\n", trueWidth)
pct <- 100 *sd(bootWidth) / trueWidth
cat("\nCV of widths in percent", pct)
The doc for runPairsBoot.R states:
"xConc, xFlux are vectors of length iboot, of the change in flow normalized concentration of flux computed by each bootstrap replicate (mg/L)"
and xConc, xFlux are computed on line 161-162 of runPairsBoot.R as:
xConc_here <- (2 * regDeltaConc) - (c22 - c11)
xFlux_here <- (2 * regDeltaFlux) - (f22 - f11)
Where regDeltaConc is the change in flow normalized flux computed by runPairs.R:
(c22 - c11)
My question is, if both xConc and regDeltaConc represent the change in flow normalized concentration, why is one: (2 * regDeltaConc) - (c22 - c11)
and the other: (c22 - c11).
Should they both be c22-c11 (and the same for flux)?
I am using the EGRETci package to evaluate water quality trends for data in the Mississippi, Minnesota, and St. Croix Rivers. We investigate trends on a calendar year basis. I was just looking for some clarification to make sure I am understanding the bootstrap results correctly.
Because I am investigating trends on a calendar year basis, I made sure for my eList that paStart=1 and paLong=12. For this site the first water quality sample was 2/8/1985 and the last sample was 12/19/2014. I then wanted to see the trend for calendar year 2005 to calendar year 2014 (last 10 years). Here are the first few lines of my bootstrap output text file:
Minnesota River at Jordan Chloride
Calendar Year
Bootstrap process, for change from Water Year 2005 to Water Year 2014
data set runs from WaterYear 1985 to Water Year 2015
Bootstrap block length in days 200
bootBreak is 39 confStop is 0.7
I am confused by how the program summarizes data by water year. The program knows I want to evaluate by calendar year, but the summary still says I am looking at water year 2005 to 2014. It then says I have samples from water year 1985 to water year 2015, because I had samples after 10/1/2014. So does the output mean I am looking at change from calendar year 2005 to 2014 but says “water year” anyways? Or did it change my request to look at water year 2005 to water year 2014 after all? I just want to make sure my 2014 October –December samples are being used and from this summary I’m not sure they are.
I took a look at this question and I think I know what is going on. My question to you is ask you to run look at your original WRTDS results for calendar year for the period of interest to you.
That means you want to do this:
eList <- setPA(paStart=1, paLong=12)
tableChange(eList, yearPoints=c(2005,2014))
from that output take note of the change values for concentration in mg/L and for flux in 10^6 kg/yr.
Now use EGRETci and do this:
caseSetUp <- trendSetUp(eList)
then it will ask you for first water year, tell it 2005, then last water year, tell it 2014.
then it will ask for nBoot, bootBreak and blockLength, You can give it 5 and 5 and 200. (It won't really matter here, I just want you to be able to get the result quickly - when you really run it you will want bigger values for the first two of these).
Then.
eBoot <- wBT(eList,caseSetUp)
when you look at the results, near the top you will have two lines that say something like this:
WRTDS estimated concentration change is xxx mg/L
WRTDS estimated flux change is yyy 10^6 kg/yr
The xxx should match the concentration change and the yyy the flux change you got in tableChange.
I tried it with an example of my own and it worked fine. Assuming the numbers work out for you then the interpretation is this. This is really a trend test on change from calendar year 2005 to calendar year 2014.
The upshot is that we need to change what our automatic output says, but the computations are doing exactly what you want. When we started designing the software we made it work just for Water Year. Then we changed to making it more flexible, but we never went back and changed the message that prints out.
Let us know if this experiment works.
Laura: can you look at a way of changing the code to reflect what is really being tested? We also needs to see that it works for any arbitrary choice of paStart and paLong so it might not even be a whole year. The paLong and paStart in eList$INFO are driving the test, but the stuff we print doesn't really reflect that. Probably need to say what the PA actually is and defines the starting and ending of the period being evaluated.
Error in hist.default(reps, breaks = xSeq, yaxs = "i", xaxs = "i", tcl = 0.5, :
some 'x' not counted; maybe 'breaks' do not span range of 'x
Here is the list:
rep Concentration | Flux
value nPos post_p Lower Upper | value nPos post_p Lower Upper
1 -0.0309 0 0.25 -0.0309 -0.0309 | 0.283 1 0.75 0.283 0.283
2 -0.0168 0 0.167 -0.0309 -0.0168 | 0.6621 2 0.833 0.283 0.6621
3 -0.03 0 0.125 -0.0309 -0.0168 | 0.3014 3 0.875 0.283 0.6621
4 -0.0181 0 0.1 -0.0309 -0.0168 | 0.7496 4 0.9 0.283 0.7496
5 -0.0169 0 0.0833 -0.0309 -0.0168 | 0.5048 5 0.917 0.283 0.7496
6 -0.0101 0 0.0714 -0.0309 -0.0101 | 0.7315 6 0.929 0.283 0.7496
7 -0.0636 0 0.0625 -0.0636 -0.0101 | -0.4792 6 0.812 -0.4792 0.7496
8 -0.0528 0 0.0556 -0.0636 -0.0101 | -0.6414 6 0.722 -0.6414 0.7496
9 -0.0171 0 0.05 -0.0636 -0.0101 | 0.6196 7 0.75 -0.6414 0.7496
10 -0.0307 0 0.0455 -0.0636 -0.0101 | 0.5275 8 0.773 -0.6414 0.7496
11 -0.013 0 0.0417 -0.0636 -0.0101 | 0.8981 9 0.792 -0.6414 0.8981
12 -0.017 0 0.0385 -0.0636 -0.0101 | 0.7489 10 0.808 -0.6414 0.8981
13 -0.0452 0 0.0357 -0.0636 -0.0101 | 0.4748 11 0.821 -0.6414 0.8981
14 -0.0345 0 0.0333 -0.0636 -0.0101 | 0.2384 12 0.833 -0.6414 0.8981
15 -0.0337 0 0.0312 -0.0636 -0.0101 | 0.2557 13 0.844 -0.6414 0.8981
16 -0.0265 0 0.0294 -0.0636 -0.0101 | 0.3773 14 0.853 -0.6414 0.8981
17 -0.0541 0 0.0278 -0.0636 -0.0101 | -0.3951 14 0.806 -0.6414 0.8981
18 -0.0287 0 0.0263 -0.0636 -0.0101 | 0.6072 15 0.816 -0.6414 0.8981
19 -0.0486 0 0.025 -0.0636 -0.0101 | -0.5431 15 0.775 -0.6414 0.8981
20 -0.0262 0 0.0238 -0.0632 -0.0103 | 0.4925 16 0.786 -0.6364 0.8907
21 -0.0859 0 0.0227 -0.0836 -0.0104 | -0.7203 16 0.75 -0.7124 0.8833
22 -0.0391 0 0.0217 -0.0825 -0.0105 | -0.0548 16 0.717 -0.7085 0.8759
23 -0.00791 0 0.0208 -0.0814 -0.00835 | 0.9121 17 0.729 -0.7045 0.9093
24 -0.042 0 0.02 -0.0803 -0.00846 | 0.878 18 0.74 -0.7006 0.9086
25 -0.047 0 0.0192 -0.0792 -0.00857 | -0.02709 18 0.712 -0.6966 0.9079
26 -0.0191 0 0.0185 -0.0781 -0.00868 | 0.5373 19 0.722 -0.6927 0.9072
27 -0.0403 0 0.0179 -0.077 -0.00879 | 0.4758 20 0.732 -0.6887 0.9065
28 -0.0671 0 0.0172 -0.0774 -0.0089 | 0.0856 21 0.741 -0.6848 0.9058
29 -0.021 0 0.0167 -0.0765 -0.00901 | 0.2997 22 0.75 -0.6808 0.9051
30 -0.0259 0 0.0161 -0.0756 -0.00912 | 0.007642 23 0.758 -0.6769 0.9044
31 -0.00347 0 0.0156 -0.0746 -0.00613 | 0.6906 24 0.766 -0.6729 0.9037
32 -0.0465 0 0.0152 -0.0737 -0.00635 | -0.6041 24 0.742 -0.669 0.903
33 -0.0117 0 0.0147 -0.0727 -0.00658 | 0.06099 25 0.75 -0.665 0.9023
34 -0.0362 0 0.0143 -0.0718 -0.0068 | 0.8749 26 0.757 -0.6611 0.9016
35 -0.00482 0 0.0139 -0.0709 -0.00455 | 0.4456 27 0.764 -0.6571 0.9009
36 -0.0547 0 0.0135 -0.0699 -0.00462 | -0.8404 27 0.743 -0.7383 0.9002
37 -0.0518 0 0.0132 -0.069 -0.00469 | -0.2126 27 0.724 -0.7323 0.8995
38 -0.0171 0 0.0128 -0.0681 -0.00476 | 0.3157 28 0.731 -0.7263 0.8988
39 -0.0405 0 0.0125 -0.0671 -0.00482 | 0.3362 29 0.738 -0.7203 0.8981
The master
branch of this repository will soon be renamed from master
to main
, as part of a coordinated change across the USGS-R and USGS-VIZLAB organizations. This is part of a broader effort across the git community to use more inclusive language. For instance, git, GitHub, and GitLab have all changed or are in the process of changing their default branch name.
We will make this change early in the week of February 6, 2022. The purpose of this issue is to give notification of the change and provide information on how to make it go smoothly.
When this change is made, forks and clones of the repository will need to be updated as well. Here are instructions for making the necessary changes.
master
to main
(<your repository> -> Settings -> Branches).First, update your local git configuration so that new repositories created locally will have the correct default branch: git config --global init.defaultBranch main
.
Now, for any forks, do the following:
master
to main
.git branch -m master main
git fetch origin
git branch -u origin/main main
git remote set-head origin -a
Consider doing it like we did with WRTDS:
https://github.com/DOI-USGS/EGRET/blob/main/R/runSurvReg.R#L67
After a call to blankTime to filter a time span with missing data, it would seem that a call to ciBands (after running the required commands prior to calling ciBands) doesn't work? The error message is shown below, it is one I haven't seen before and I think it may be related to the fact that I've applied blankTime higher up in the script:
...
CIAnnualResults <- ciBands(eList, repAnnual, probs)
#Error in quantile.default(manyAnnualResults[iYear, 1, 1:nBoot], prob = probs, :
# missing values and NaN's not allowed if 'na.rm' is FALSE
Specifically, define the stuff coming out of bootOut
For calculating the right slice of "surface"s
The function wBT()
returns a list object (eBoot). The vectors "pFlux" and "pConc" are the reported change in flow-normalized values expressed in percent change. However, the documentation that describes what this functions returns does not properly render the percent symbol.
https://github.com/USGS-R/EGRETci/blob/d1c4845e53a8599faef68015596ce20dc5cd316d/R/wBTCode.R#L226
library(EGRETci)
?wBT
If you ask for wBT
for years past the model results, you get a cryptic message.
It takes about 10 minutes for wBT to run and I haven't successfully gotten ciCalculations to run. Any ideas?
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.