116 lines
7.1 KiB
R
116 lines
7.1 KiB
R
library(tidyverse)
|
|
library(fixest)
|
|
library("geosphere")
|
|
|
|
DF <- readRDS("Data/Output_Data/Full_Data_Set.rds") %>% mutate(year=as.numeric(year)) %>% group_by(wdid) %>% mutate(SBD1=max(SBD1)) %>% ungroup
|
|
TEMP <- read_csv("Data/Output_Data/Div3_Pumping_Data.csv")
|
|
TEMP %>% group_by(wdid,year) %>% filter(n()>1) %>% print(n=100)
|
|
#EMPTY_START <- DF %>% filter(year<=2010) %>% group_by(wdid) %>% filter(sum(AF)==0) %>% select(wdid) %>% unique
|
|
#DF <- DF %>% anti_join(EMPTY_START)
|
|
OLD <- readRDS("Test/df_half.Rds") %>% mutate(wdid=GW_wdid,SBD1_OLD=SBD1,year=Year,AF_OLD=AF) %>% select(wdid,year,SBD1_OLD,AF_OLD)
|
|
BOTH <- DF %>% select(wdid,year,SBD1,AF) %>% full_join(OLD) %>% filter(year<2020)
|
|
BOTH %>% filter(is.na(SBD1))
|
|
BOTH %>% filter(AF!=AF_OLD) %>% print(n=800)
|
|
BOTH %>% filter(AF!=AF_OLD)
|
|
OLD %>% filter(wdid=='2405396')
|
|
BOTH %>% filter(wdid=='2405396')
|
|
|
|
BOTH %>% filter(is.na(SBD1_OLD))
|
|
|
|
|
|
|
|
DF$SBD1 <- ifelse(DF$SBD1_year!=Inf,1,0)
|
|
#DF <- DF %>% select(-ditch_2200627,-ditch_2200541,-ditch_3500570)
|
|
DF$SBD1_year <- ifelse(DF$SBD1_year==2009,2011,DF$SBD1_year)
|
|
DF$POST <- ifelse(DF$year>=DF$SBD1_year,1,0)
|
|
DF <- DF %>% left_join(DF %>% group_by(wdid,CREP_any) %>% summarize(CREP_year=min(year)) %>% filter(CREP_any==1) %>% select(-CREP_any)) %>% mutate(CREP_year=ifelse(is.na(CREP_year),Inf,CREP_year))
|
|
#feols(AF~SBD1*POST+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028|wdid+SBD2+SBD3+SBD4+SBD5+SBD6+wdid+ditch_2000812^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year+ditch_2200627^year+ditch_2200541^year+ditch_3500570^year,data=DF)
|
|
###########
|
|
#DF$SBD1_year <- ifelse(DF$SBD1_year==Inf,10000,DF$SBD1_year)
|
|
|
|
|
|
|
|
DF$POST_TEST <- ifelse(DF$year>=2011,1,0)
|
|
DF %>% filter(year==2009) %>% pull(wdid) %>% unique %>% length
|
|
TEST_DATA <- DF %>% mutate(SBD1=ifelse(SBD1==1 & SBD1_year<2020,1,0))
|
|
SUN_MOD <- feols(AF~SBD1+sunab(SBD1_year,year)+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028+per_alfalfa+per_potatoes|wdid+year+SBD2^year+SBD3^year+SBD4^year+SBD5^year+SBD6^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year,data=DF %>% filter(SBD1_year==Inf|SBD1_year==2011) )
|
|
etable(feols(AF~SBD1*POST_TEST|wdid+year,data=TEST_DATA %>% filter(year<2019)))
|
|
|
|
coefplot(SUN_MOD)
|
|
SUN_MOD <- feols(AF~SBD1+sunab(SBD1_year,year)+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028+per_alfalfa+per_potatoes|wdid+year+SBD2^year+SBD3^year+SBD4^year+SBD5^year+SBD6^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year,data=DF )
|
|
|
|
TEST <- DF %>% group_by(year,contacts,SBD1) %>% summarize(AF=sum(AF),SBD1_year=min(SBD1_year),POST=max(POST),CREP_temp=max(CREP_temp),CREP_perm=max(CREP_perm),SBD1_temp_fallow=max(SBD1_temp_fallow),SB22_028=max(SB22_028),SBD_half_fallow=max(SBD_half_fallow),SBD1_purchase=max(SBD1_purchase),SBD2=max(SBD2),SBD3=max(SBD3),SBD4=max(SBD4),SBD5=max(SBD5),SBD6=max(SBD6),ditch_2000631=max(ditch_2000631),ditch_Other=max(ditch_Other),ditch_2000829=max(ditch_2000829),ditch_2000798=max(ditch_2000798),ditch_2000816=max(ditch_2000816),ditch_2000753=max(ditch_2000753),ditch_2000753=max(ditch_2000753),ditch_2000623=max(ditch_2000623),CREP_any=max(CREP_any)) %>% ungroup
|
|
TEST %>% arrange(contacts,year)
|
|
TEST
|
|
TEST %>% select(contacts,year,SBD1,POST) %>% arrange(contacts,year) %>% filter(SBD1==1)
|
|
TEST_REG <- feols(log(AF+0.001)~sunab(SBD1_year,year)+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028|contacts+SBD2^year+SBD3^year+SBD4^year+SBD5^year+SBD6^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year,TEST )
|
|
TEST_REG <- feols(AF~sunab(SBD1_year,year)+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028|contacts,TEST )
|
|
|
|
coefplot(TEST_REG)
|
|
TEST
|
|
etable(TEST_REG ,agg="cohort")
|
|
|
|
coefplot(SUN_MOD)
|
|
#coefplot(feols(AF~sunab(CREP_year,year)|wdid+year,data=DF))
|
|
coefplot(feols(AF~sunab(CREP_year,year)|wdid+year+ditch_2000812^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year+ditch_2200627^year+ditch_2200541^year+ditch_3500570^year,data=DF))
|
|
|
|
|
|
|
|
|
|
|
|
|
|
GET_PROGRAM <-function(PROGRAM_NUMBER,DATA,CUTOFF=0.25,YEAR_RANGE=2009:2025){
|
|
#DATA <- DF
|
|
#PROGRAM_NUMBER <- 4
|
|
#YEAR_RANGE=2009:2025
|
|
#CUTOFF <- 0.25
|
|
DATA_ORIG <- DATA
|
|
MASTER_WELL_LIST <- DATA %>% select(wdid,longitude,latitude) %>% unique
|
|
WELLS_IN_PROGRAM <- DATA[DATA[,PROGRAM_NUMBER]==1,]
|
|
PROGRAM_START <- min(WELLS_IN_PROGRAM$year)
|
|
PROGRAM_YEAR_RANGE <- min(WELLS_IN_PROGRAM$year):max(YEAR_RANGE)
|
|
for(YEAR in PROGRAM_YEAR_RANGE){
|
|
# YEAR <-2014
|
|
C_PROGRAM_WELLS <- WELLS_IN_PROGRAM[WELLS_IN_PROGRAM[,"year"]==YEAR,c("longitude","latitude")] %>% unique %>% as.matrix
|
|
RES_DISTANCE <- c()
|
|
for(i in 1:nrow(MASTER_WELL_LIST)){
|
|
RES_DISTANCE[i] <-ifelse((2.38084*min(distHaversine(MASTER_WELL_LIST[i,c("longitude","latitude")],C_PROGRAM_WELLS))/5280)<=CUTOFF,1,0)
|
|
}
|
|
C_RES <- cbind(MASTER_WELL_LIST[,1],RES_DISTANCE) %>% as_tibble
|
|
colnames(C_RES) <-c("wdid",paste0("close_",colnames(DATA_ORIG)[PROGRAM_NUMBER]))
|
|
C_RES$year <- YEAR
|
|
if(exists("RES")){RES <- rbind(RES,C_RES)} else{RES <- C_RES}
|
|
}
|
|
if(PROGRAM_START>min(YEAR_RANGE)){
|
|
FILL_IN_YEARS <- min(YEAR_RANGE):(PROGRAM_START-1)
|
|
length(FILL_IN_YEARS )
|
|
5*(nrow(MASTER_WELL_LIST))
|
|
|
|
MASTER_WELL_LIST[,1]
|
|
FILL_IN_WELLS <- rep(t(MASTER_WELL_LIST[,1]),length(FILL_IN_YEARS))
|
|
ZEROS <- rep(0,length(FILL_IN_WELLS))
|
|
FILL_IN_YEARS <- rep( FILL_IN_YEARS,nrow(MASTER_WELL_LIST)) %>% sort
|
|
FILL_IN <- cbind(FILL_IN_WELLS,ZEROS,FILL_IN_YEARS) %>% as_tibble
|
|
FILL_IN <- FILL_IN %>% mutate(ZEROS=as.numeric(ZEROS),FILL_IN_YEARS=as.numeric(FILL_IN_YEARS))
|
|
|
|
colnames(FILL_IN) <- c("wdid",paste0("close_",colnames(DATA_ORIG)[PROGRAM_NUMBER]),"year")
|
|
|
|
RES <- rbind(FILL_IN,RES) %>% unique
|
|
}
|
|
RES <- RES[,c(1,3,2)]
|
|
return(RES)
|
|
}
|
|
|
|
DF <- DF %>% left_join(GET_PROGRAM(4,DF)) %>%left_join(GET_PROGRAM(5,DF)) %>% left_join(GET_PROGRAM(6,DF) )
|
|
CLOSE_YEAR <-
|
|
CLOSE_YEAR <-DF %>% filter(close_CREP_any==1) %>% group_by(wdid) %>% summarize(close_CREP_year=min(year)) %>% unique
|
|
DF <- DF %>% left_join(CLOSE_YEAR) %>% mutate(close_CREP_year=ifelse(is.na(close_CREP_year),Inf,close_CREP_year))
|
|
|
|
|
|
summary(DF$close_CREP_year)
|
|
coefplot(feols(AF~sunab(close_CREP_year,year)+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028|wdid+SBD1^year+SBD2^year+SBD3^year+SBD4^year+SBD5^year+SBD6^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year,data=DF ))
|
|
|
|
etable(feols(AF~sunab(close_CREP_year,year)+CREP_temp+CREP_perm+SBD_half_fallow+SBD1_purchase+SBD1_temp_fallow+SB22_028|wdid+SBD1^year+SBD2^year+SBD3^year+SBD4^year+SBD5^year+SBD6^year+ditch_2000631^year+ditch_Other^year+ditch_2000829^year+ditch_2000798^year+ditch_2000816^year+ditch_2000753^year+ditch_2000753^year+ditch_2000623^year,data=DF ),agg="cohort")
|
|
|
|
|