Skip to content

Instantly share code, notes, and snippets.

View SwampThingPaul's full-sized avatar

Paul Julian SwampThingPaul

View GitHub Profile
# This script relies on the magick library, see this vignette for more information:
# https://cran.r-project.org/web/packages/magick/vignettes/intro.html
#
# If you have Windows or Mac OS, I believe ImageMagick STL is integrated with the CRAN distribution
# of the R magick library. If on Linux, see the "Build from source" section of the above URL
#
# This script will add year annotations from GIF animations generated and downloaded from the
# LT-GEE Time Series Animator open access Earth Engine App:
# https://emaprlab.users.earthengine.app/view/lt-gee-time-series-animator
#
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
plt.style.use('ggplot')
fig, (ax1, ax2) = plt.subplots(ncols=2,figsize=(16,6))
plt.xlim((0, 10))
plt.ylim((0, 7))
plt.tight_layout(w_pad=1.5)
#red line
N=60
mu=0
sd=2
np.random.seed(0)
ran = np.random.normal(size=N)
error1 = sd**2 * ran + mu
error2 = sd*.5 * ran + mu
lin = np.linspace(-15., 15., num=N)
import scipy.odr as odr
def odr_line(B, x):
y = B[0]*x + B[1]*x**2
return y
def perform_odr(x, y, xerr, yerr):
quadr = odr.Model(odr_line)
mydata = odr.Data(x, y, wd=1./xerr, we=1./yerr)
#mydata = odr.Data(x, y)
import numpy.linalg as la
def tls(X,y):
if X.ndim is 1:
n = 1 # the number of variable of X
X = X.reshape(len(X),1)
else:
n = np.array(X).shape[1]
Z = np.vstack((X.T,y)).T
import numpy as np
from numpy.linalg import inv
import statsmodels.api as sm
# from scratch
x = sm.add_constant(x) # add constant in the 0 index
b = inv(x.T.dot(x)).dot(x.T).dot(y)
yest_ols = np.array([b[2]*v**2 + b[1]*v + b[0] for v in x.T[0]])
# with using numpy.linalg.lstsq
library(jsonlite);library(sf);library(sp);library(geojsonio)
watershed = function(state, lon, lat, sf=TRUE){
p1 = 'https://streamstats.usgs.gov/streamstatsservices/watershed.geojson?rcode='
p2 = '&xlocation='
p3 = '&ylocation='
p4 = '&crs=4326&includeparameters=false&includeflowtypes=false&includefeatures=true&simplify=true'
query <- paste0(p1, state, p2, toString(lon), p3, toString(lat), p4)
mydata <- fromJSON(query, simplifyVector = FALSE, simplifyDataFrame = FALSE)
poly_geojsonsting <- toJSON(mydata$featurecollection[[2]]$feature, auto_unbox = TRUE)
@menugget
menugget / plot.stacked.2.R
Last active February 25, 2019 18:38
Stacked plot
#plot.stacked makes a stacked plot where each y series is plotted on top
#of the each other using filled polygons
#
#Arguments include:
#'x' - a vector of values
#'y' - a matrix of data series (columns) corresponding to x
#'order.method' = c("as.is", "max", "first")
# "as.is" - plot in order of y column
# "max" - plot in order of when each y series reaches maximum value
# "first" - plot in order of when each y series first value > 0
@menugget
menugget / plot.stream.4.R
Last active February 25, 2019 18:46
Stream plot.
#plot.stream makes a "stream plot" where each y series is plotted
#as stacked filled polygons on alternating sides of a baseline.
#
#Arguments include:
#'x' - a vector of values
#'y' - a matrix of data series (columns) corresponding to x
#'order.method' = c("as.is", "max", "first")
# "as.is" - plot in order of y column
# "max" - plot in order of when each y series reaches maximum value
# "first" - plot in order of when each y series first value > 0
set.seed(2)
x <- 1:100
y <- 20 + 3 * x
e <- rnorm(100, 0, 60)
y <- 20 + 3 * x + e
plot(x,y)
yx.lm <- lm(y ~ x)
lines(x, predict(yx.lm), col="red")