Meteorological Review

Show the code

#Tidyverse
library(plyr); library(dplyr)
library(tidyverse); library(lubridate);library(magrittr)
#Time series
library(zoo);library(xts);library(timetk);library(hydroTSM);
library(dygraphs)
#ggplot support
library(ggmap);library(ggpmisc);library(ggrepel);library(ggsn)
#spreadsheets
library(openxlsx);library(readxl)
#paralel computing
library(parallel);library(pbapply)
#general support
library(here);library(glue);library(udunits2);library(usethis);
library(matrixStats);library(pander);library(reshape2);library(Hydro)
library(geosphere);library(nasapower)
#data table support
library(data.table);library(dtplyr);library(arrow)
#raster support
library(raster);library(sp);library(kgc)


#These directories are implemented to organize the regional analysis
use_directory("data")
use_directory("data/csv")
use_directory("data/netcdf")
use_directory("data/rds")
use_directory("graphs")
use_directory("reports")
use_directory("scripts")
use_directory("outputs")

glue::glue("Work directory: {here()}")

Site Location

Show the code

Longitude= -162.85
Latitude = 68.07

source(here::here("scripts","hydro_support.r"))
#For more information


#Get elevation from google maps
Elevation=GM_elev(Longitude,Latitude,key_google = key_google)

site.info<-data.frame(Location="Site",
                      Longitude=Longitude,
                      Latitude=Latitude,
                      Elevation=Elevation) 

#Site in looped with variables
site.info_xy<-site.info %>% 
  dplyr::rename(x="Longitude",y="Latitude",z="Elevation") %>% 
  reshape2::melt(id.vars="Location")

pandoc.table(site.info)

Location	Longitude	Latitude	Elevation
Site	-162.8	68.07	299.7

Show the code

#minimum year of information used, to be defined after seen the available info
min.ywi<-10 

Downloading records

Show the code

source(here::here("scripts","usgs_download.r"))

USGS_flow_index(Latitude = Latitude,Longitude = Longitude,n = 30)

USGS_flow_download(x = NAflow.idx$id,monthly.zero = T) %>% invisible()

[1] “Stations N 1 : 15746983 / RED DOG MINE CLEAN WATER D NR KIVALINA AK” [1] “Stations N 2 : 15746988 / NF RED DOG C NR KIVALINA AK” [1] “Stations N 3 : 15746990 / RED DOG C AB MOUTH NR KIVALINA AK” [1] “Stations N 4 : 15746980 / IKALUKROK C AB RED DOG C NR KIVALINA AK” [1] “Stations N 5 : 1574699020 / IKALUKROK C 0.6 MI BL RED DOG C NR KIVALINA AK” [1] “Stations N 6 : 15746991 / IKALUKROK C BL RED DOG C NR KIVALINA AK” [1] “Stations N 7 : 15746900 / WULIK R AB FERRIC C NR KIVALINA AK” [1] “Stations N 8 : 15747000 / WULIK R BL TUTAK C NR KIVALINA AK” [1] “Stations N 9 : 15746000 / NOATAK R AT NOATAK AK” [1] “Stations N 10 : 15748000 / OGOTORUK C NR POINT HOPE AK” [1] “Stations N 11 : 15743000 / JUNE C NR KOTZEBUE AK” [1] “Stations N 12 : 15744500 / KOBUK R NR KIANA AK” [1] “Stations N 13 : 15744000 / KOBUK R AT AMBLER AK” [1] “Stations N 14 : 15716010 / HUMBOLT C NR SERPENTINE HOT SPRINGS NR NOME AK” [1] “Stations N 15 : 15743850 / DAHL C NR KOBUK AK” [1] “Stations N 16 : 15712000 / KUZITRIN R NR NOME AK” [1] “Stations N 17 : 15803000 / MEADE R AT ATKASUK AK” [1] “Stations N 18 : 15668200 / CRATER C NR NOME AK” [1] “Stations N 19 : 15635000 / ELDORADO C NR TELLER AK” [1] “Stations N 20 : 15583500 / ETTA C NR COUNCIL AK” [1] “Stations N 21 : 15625850 / STEWART R 0.1 MI BL BOULDER C MTH NR NOME AK” [1] “Stations N 22 : 15625900 / STEWART R 0.2 MI BL DURRANT C MTH NR NOME AK” [1] “Stations N 23 : 15621000 / SNAKE R NR NOME AK” [1] “Stations N 24 : 15798700 / NUNAVAK C NR BARROW AK” [1] “Stations N 25 : 15799000 / ESATKUAT C NR BARROW AK” [1] “Stations N 26 : 15799300 / ESATKUAT LAGOON OUTLET AT BARROW AK” [1] “Stations N 27 : 15564900 / KOYUKUK R AT HUGHES AK” [1] “Stations N 28 : 15830000 / MIGUAKIAK R NR TESHEKPUK LK NR LONELY AK” [1] “Stations N 29 : 15565200 / YUKON R NR KALTAG AK” [1] “Stations N 30 : 15875000 / COLVILLE R AT UMIAT AK”

Show the code

#selected the station with more than 20 values
gauge_sel<-lapply(NAflow.db$annual.peak,FUN=function(x)length(x))>20

flow_peak_tbl<-lapply(NAflow.db$annual.peak[gauge_sel],FUN=function(x){fortify.zoo(x)}) %>% 
  do.call(rbind,.) %>% 
  rownames_to_column(var="id") %>% 
  mutate(id=str_split(id,"\\.",simplify = T) %>% .[,1] %>% str_remove_all("X")) %>% 
  mutate(mon_day=month(Index)+day(Index)/31,
         yr=year(Index)) %>% 
  rename(date="Index",peak="x")

flow_peak_idx<-readNWISsite(siteNumbers =unique(flow_peak_tbl$id)) %>%
  dplyr::select(id="site_no",name="station_nm",
                lat="dec_lat_va",lon="dec_long_va",area_mi2="drain_area_va") %>% 
  mutate(area_km2=ud.convert(area_mi2,"mi2","km2")) %>% 
  mutate(dist_km=geosphere::distVincentyEllipsoid(p1=cbind(lon,lat),
                                               p2=cbind(Longitude,Latitude))/1000) %>% 
  arrange(dist_km)

Gauge results

Show the code

pandoc.table(flow_peak_idx,
             split.table=Inf,round=2)

id	name	lat	lon	area_mi2	area_km2	dist_km
15746988	NF RED DOG C NR KIVALINA AK	68.08	-162.9	14	36.26	2.07
15747000	WULIK R BL TUTAK C NR KIVALINA AK	67.88	-163.7	705	1826	40.88
15744500	KOBUK R NR KIANA AK	66.97	-160.1	9480	24553	168.5
15743850	DAHL C NR KOBUK AK	66.95	-156.9	10.9	28.23	282.8
15668200	CRATER C NR NOME AK	64.93	-164.9	22	56.98	361.6
15635000	ELDORADO C NR TELLER AK	64.96	-166.2	6	15.54	377.6
15621000	SNAKE R NR NOME AK	64.56	-165.5	86	222.7	408.7
15798700	NUNAVAK C NR BARROW AK	71.26	-156.8	2.9	7.51	426.1
15564900	KOYUKUK R AT HUGHES AK	66.05	-154.3	17990	46594	436
15875000	COLVILLE R AT UMIAT AK	69.36	-152.1	13860	35897	457.3

Figures

Peak flow distribution

Show the code

flow_peak_idx<-flow_peak_idx %>% 
  mutate(peak_type=if_else(id%in%c("15746988","15747000","15635000","15743850"),
                           "rain","snow")) %>% 
  mutate(peak_type=if_else(id%in%c("15668200","15621000","15564900"),"mixed",peak_type))

left_join(flow_peak_tbl,flow_peak_idx) %>% 
  mutate(label=paste(name,"- area:",round(area_km2,0),"km2") %>% 
           str_wrap(20)) %>% 
  dplyr::filter(mon_day>3) %>% 
  ggplot(data=.,aes(x=mon_day,y=peak,col=yr,shape=as.factor(peak_type)))+
  geom_point(size=2.5,alpha=0.8)+
  scale_color_binned(type = "viridis")+
  theme_light()+
  scale_x_continuous(breaks=1:12,labels=month.abb)+
  labs(x="peak month",y="Peak Magnitude [m3/s]",shape="peak type:",
       title="Month location and peak flow magnitude in a calendar year")+
  facet_wrap(~label,scales = "free_y")+
  theme(strip.text.x = element_text(size=8),
        plot.title=element_text(size=16))+
  theme(legend.position="bottom")

Show the code

left_join(flow_peak_tbl,flow_peak_idx) %>% 
  mutate(label=paste(name,"- area:",round(area_km2,0),"km2") %>% 
           str_wrap(20)) %>%
  mutate(mon=month(date)) %>% 
  dplyr::filter(mon_day>3) %>% 
  ggplot(data=.,aes(x=mon,fill=as.factor(peak_type)))+
  geom_histogram()+
  facet_wrap(~label)+
  labs(fill="peak type",title="Frequency of peak flows by month")+
  scale_x_continuous(breaks=1:12,labels=month.abb)+
  theme_light()+
    theme(strip.text.x = element_text(size=8),
        plot.title=element_text(size=16))+
  theme(legend.position="bottom")

Overall map

Show the code

#used for ggsn graphical scale, same ranges but 8% smaller
xy <-ggmap::make_bbox(lon = c(flow_peak_idx$lon,Longitude),
                      lat = c(flow_peak_idx$lat,Latitude),f=0.1) 

#ggmap capture for goggle maps
map <- get_map(xy,scale=2,zoom=5,maptype="hybrid",source="google",
               messaging=FALSE)

ggmap(map)+
  geom_point(data=site.info,
             aes(x=Longitude,y=Latitude),shape=5,col="red",size=5)+
  geom_label_repel(data=site.info,
                   aes(x=Longitude,y=Latitude,label=Location))+
  geom_point(data=flow_peak_idx,aes(x=lon,y=lat,size=area_km2,col=as.factor(peak_type)))+
  geom_label_repel(data=flow_peak_idx,aes(x=lon,y=lat,label=name %>% 
                                            str_wrap(12)),size=3.5)+
  labs(col="peak type:",size="area km2:")