---
title: "Trajectories of Brownian motion"
output: html_document
---

```{r}
rm(list=ls())
library(berryFunctions)
library(ggplot2)
library("gridExtra")

sim.bm<-function(x,n_path=n_path){
  tau<-x[1]
  bm<-cumsum(c(0,tau*rnorm(n_path-1)))
  return(bm)
}

##### tau=1
set.seed(463)
n_path<-100
start<-1:(n_path-1)
end<-2: n_path
reps<-100

tau1<-1
x<-c(tau1)
sim.bm.array1<-array(NA,c(reps,n_path))
for(i in 1:reps){
  sim.bm.array1[i,]<-sim.bm(x,n_path=n_path)
}

tau2<-4
x<-c(tau2)
sim.bm.array2<-array(NA,c(reps,n_path))
for(i in 1:reps){
  sim.bm.array2[i,]<-sim.bm(x,n_path=n_path)
}

p1<-ggplot() + ylim(range(unlist(sim.bm.array1),unlist(sim.bm.array2)))+xlab("Time Steps") + ylab("Trait Values")
p2<-ggplot() + ylim(range(unlist(sim.bm.array1),unlist(sim.bm.array2)))+xlab("Time Steps") + ylab("Trait Values")

for(repIndex in 1:reps){
  #repIndex<-1
  df1 <- data.frame(x = start, y = c(sim.bm.array1[repIndex,start]), xend=end, yend=c(sim.bm.array1[repIndex,end])) 
  p1<- p1 +  geom_segment(data=df1, mapping=aes(x=x, y=y, xend=xend, yend=yend), size=0.5, color="salmon")
}

dim(sim.bm.array1)
df11<-data.frame(x=rep(reps,n_path),y=sim.bm.array1[,n_path])
p1<- p1 +  geom_point(data=df11, mapping=aes(x=x, y=y), size=1, color="deepskyblue")
#p1

# sim.bm.ave1<-apply(sim.bm.array1,2,mean)
# df.ave1 <- data.frame(x = start, y = c(sim.bm.ave1[start]), xend=end, yend=c(sim.bm.ave1[end])) 

#p1<- p1 +  geom_segment(data=df.ave1, mapping=aes(x=x, y=y, xend=xend, yend=yend), size=1, color="deepskyblue")

p1<-p1+scale_color_manual(values=c("#CC6666"))
p1<-p1+theme(axis.title.x=element_blank(),
             axis.text.x=element_blank(),
             axis.ticks.x=element_blank(),
             axis.title.y=element_blank()
)# + ggtitle("BM")
p1=p1+theme_bw()+coord_flip()



hist_p1 <- ggplot()+theme_bw()+geom_histogram(aes(  sim.bm.array1[,n_path]  ),color="navy", fill="deepskyblue")+ xlim(range(unlist(sim.bm.array1[,n_path]),unlist(sim.bm.array2[,n_path]))) + theme(axis.title.x=element_blank(),axis.text.x=element_blank(),axis.ticks.x=element_blank(),axis.text.y=element_blank(),axis.ticks.y=element_blank(),axis.title.y=element_blank())
#hist_p1
##### tau=4

for(repIndex in 1:reps){
  df2 <- data.frame(x = start, y = c(sim.bm.array2[repIndex,start]), xend=end, yend=c(sim.bm.array2[repIndex,end])) 
  p2<- p2 +  geom_segment(data=df2, mapping=aes(x=x, y=y, xend=xend, yend=yend), size=0.5, color="salmon")
}

df21<-data.frame(x=rep(reps,n_path),y=sim.bm.array2[,n_path])
p2<- p2 +  geom_point(data=df21, mapping=aes(x=x, y=y), size=1, color="deepskyblue")
#p2


# sim.bm.ave2<-apply(sim.bm.array2,2,mean)
# df.ave2 <- data.frame(x = start, y = c(sim.bm.ave2[start]), xend=end, yend=c(sim.bm.ave2[end])) 

#p2<- p2 +  geom_segment(data=df.ave2, mapping=aes(x=x, y=y, xend=xend, yend=yend), size=1, color="deepskyblue")

p2<-p2+scale_color_manual(values=c("#CC6666"))
p2<-p2+theme(axis.title.x=element_blank(),
             axis.text.x=element_blank(),
             axis.ticks.x=element_blank(),
             axis.title.y=element_blank()
)# + ggtitle("BM")

p2=p2+theme_bw()+coord_flip()


hist_p2 <- ggplot()+theme_bw()+geom_histogram(aes(  sim.bm.array2[,n_path]  ),color="navy", fill="deepskyblue") + xlim(range(unlist(sim.bm.array1[,n_path]),unlist(sim.bm.array2[,n_path])))+ theme(axis.title.x=element_blank(),axis.text.x=element_blank(),axis.ticks.x=element_blank(),axis.text.y=element_blank(),axis.ticks.y=element_blank(),axis.title.y=element_blank())

grid.arrange(hist_p1, hist_p2, p1, p2, nrow=2, ncol = 2,widths=c(2, 2),heights=c(1,3))


```