rm(list = ls())

setwd("C:/Users/screel/Desktop/KENYA/ungulate_transects_in_R") 
# you will need to set the working directory appropriately for your machine
# recall that there is a menu option under session in RStudio

library(unmarked)

dists <- read.csv("knp.allsp.oct2012.csv", header=TRUE)
dists$trans <- as.factor(dists$trans)                       #convert transects into a categorical variable

# look at the data read into dists

dim(dists)          # dimensions (observations = rows, the 1st number; variables = columns, the 2nd number)
head(dists)         # display first 10 rows on console
tail(dists)         # display last 10 rows on console

dists.sub <- subset(dists, species == 'PU')
dim(dists.sub)
head(dists.sub)
tail(dists.sub)

dists.order <- dists.sub[order(dists.sub$dist), ]        #this orders the data from small to large distances
dists.order$dist   #display the ordered distances on the console and check dimensions again to be sure all OK
dim(dists.order)

dists.trunc <-  dists.order[which(dists.order$dist < 300), ]   
dim(dists.trunc)   #still 54 observations after truncation of the data to remove sparse information beyond 300m                                             
dists.trunc 

yDat <- formatDistData(dists.trunc, distCol="dist", transectNameCol="trans", 
                       dist.breaks=c(0, 50, 100, 150, 200, 250, 300))  #create distance categories (m). Could change these.
yDat
nrow(yDat)          #Should be  76 rows for this subset of KNP data, equal to ~120 km of transects sampled on one occasion.

covs <- read.csv("C:/Users/screel/Desktop/KENYA/ungulate_transects_in_R/knp_covs_oct2012.csv", header=TRUE)
head(covs)
tail(covs)
nrow(covs)

umf <- unmarkedFrameDS(y=as.matrix(yDat), siteCovs=covs, survey="line", dist.breaks=c(0, 50, 100, 150, 200, 250, 300), tlength=covs$length,  unitsIn="m")
head(umf)  #view the dataframe you just created, first 10 rows only.  yDat and covs are now linked.

hist(umf, xlab="distance (m)", main="", cex.lab=0.8, cex.axis=0.8)

m.half <- distsamp(~1 ~1, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.haz <- distsamp(~1 ~1, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.uni <- distsamp(~1 ~1, umf, keyfun="uniform", output="density", unitsOut="kmsq")

m.half
m.haz
m.uni


#half-norm:
m.half.1.type <- distsamp(~1 ~type, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.half.1.hab <- distsamp(~1 ~hab, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.half.1.ht <- distsamp(~1 ~ht, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.half.1.col <- distsamp(~1 ~col, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.half.1.burn <- distsamp(~1 ~burn, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.half.1.lagoon <- distsamp(~1 ~lagoon, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")
m.half.1.seaswater <- distsamp(~1 ~seas.water, umf, keyfun="halfnorm", output="density", unitsOut="kmsq")



#hazard:
m.haz.1.type <- distsamp(~1 ~type, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.haz.1.hab <- distsamp(~1 ~hab, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.haz.1.ht <- distsamp(~1 ~ht, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.haz.1.col <- distsamp(~1 ~col, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.haz.1.burn <- distsamp(~1 ~burn, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.haz.1.lagoon <- distsamp(~1 ~lagoon, umf, keyfun="hazard", output="density", unitsOut="kmsq")
m.haz.1.seaswater <- distsamp(~1 ~seas.water, umf, keyfun="hazard", output="density", unitsOut="kmsq")

fmList <- fitList(m.half=m.half, m.uni=m.uni, 
                  m.half.1.type=m.half.1.type, m.half.1.hab=m.half.1.hab, m.half.1.ht=m.half.1.ht,
                  m.half.1.col=m.half.1.col, m.half.1.burn=m.half.1.burn, m.half.1.lagoon=m.half.1.lagoon,
                  m.half.1.seaswater=m.half.1.seaswater, m.haz.1.type=m.haz.1.type, m.haz.1.hab=m.haz.1.hab,
                  m.haz.1.ht=m.haz.1.ht, m.haz.1.burn=m.haz.1.burn,
                  m.haz.1.lagoon=m.haz.1.lagoon, m.haz.1.seaswater=m.haz.1.seaswater, m.haz.1.col=m.haz.1.col, m.haz=m.haz) 
modSel(fmList)