Skip to content

Instantly share code, notes, and snippets.

View SwampThingPaul's full-sized avatar

Paul Julian SwampThingPaul

View GitHub Profile
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 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
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 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)
@SwampThingPaul
SwampThingPaul / spatial_package_install.r
Last active January 7, 2019 13:36
Package installer helper
packs=c("sp","rgdal","gstat","raster","spatstat","maptools","rgeos","tmap","GISTools","rasterVis","spdep","spsurvey")
for(i in 1:length(packs)){
test=packs[i] %in% rownames(installed.packages())
if(test==T){print("package already installed")}else{install.packages(packs[i])}
}
@SwampThingPaul
SwampThingPaul / thessian_clip.r
Last active January 7, 2019 14:37
Geospatial analysis - Thessian Polygon in R
require(spstat)
point.dat;# a shapefile of points with associated data.
study.area;# a shapefile of the study area/sampling area.
# Generate Thessian polygon and assign CRS
th=as(dirichlet(as.ppp(point.dat)),"SpatialPolygons")
proj4string(th)=proj4string(point.dat)
# Join thessian polygon with actual data
@SwampThingPaul
SwampThingPaul / InsetMap.r
Last active January 7, 2019 14:47
mapping in rstats
map2=base.map+tm_shape(TN.GM)+
tm_symbols(size=0.5,col="Geomean",breaks=c(-Inf,0.5,1,2,Inf),showNA=T,palette=cols.rmp,
title.col="Annual Geometric \nMean TN \nConcentration (mg/L)",
labels=c("\u003C 0.5","0.5 - 1.0","1.0 - 2.0", "\u003E2.0"),
border.lwd=0.5,colorNA = "white")+
tm_compass(type="arrow",position=c("left","bottom"))+
tm_scale_bar(position=c("left","bottom"))+
tm_layout(bg.color=cols[2],fontfamily = "serif",legend.outside=T,scale=1,asp=NA,
outer.margins=c(0.005,0.01,0.005,0.01),inner.margins = 0,between.margin=0,
legend.text.size=1,legend.title.size=1.25)
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
@SwampThingPaul
SwampThingPaul / HURDAT_mapping.r
Last active January 12, 2019 14:01
Hurricane Data Analysis
## Code was compiled by Paul Julian
## contact infor: pjulian@ufl.edu
#Libraries
library(HURDAT)
library(plyr)
library(sp)
library(tmap)
#Projection
@SwampThingPaul
SwampThingPaul / correlation_power.r
Created January 23, 2019 15:00
Power analysis FUN
library(pwr)
##
r.val=seq(-0.9,-0.01,0.05)
n.val=4:20
power.rslt=data.frame()
for(j in 1:length(n.val)){