# This is a script to simulate an experiment of multiple tosses of a coin,
# or, as an equivalent experiment, of one-dimensional random walk.
# The coin is tossed 'n_toss' times and this experiment is repeated 'n_rep' times.
# Heads are represented as -1 and tails as +1. Their accumulated sum is denoted by 'walk_vec'.
# Plots generated represent the one-dimensional random walk.
# Most of the time, random variable y is rather far away from its expectation of zero.
# However, the average distance ('walk_vec' divided by the number of tosses) converges to zero.
# Plots of 'walk_vec' and the average distance is displayed after every 'intermed_plots' tosses.
# Note that average distance plot is scaled up by a factor of 10.
# The distances 'walk_vec' are expected to be an arcsine distribution.
# To show this, after 'n_rep' experiments with 'n_toss' tosses, a histogram of these distances is plotted.
# The density function of an arcsine distributed random variable is plotted as a continuous line
# on the histogram to compare.

# The script uses packages VaRES and ggplot2; thus, these should be installed prior to running the script.
require(VaRES)
require(ggplot2)
require(profvis)

# Input values to be specified by the user

n_toss <- 1000             # Number of tosses in one experiment
n_rep <- 1000              # Number of repetitions of the experiment
intermed_plots <- 100      # Number of repetitions after which to display intermediate results (for no intermediate plots: 0)
                           # (criterion: intermed_plots = 0 (mod n_rep))
halt <- FALSE               # Boolean that describes whether or not to halt after every intermediate graph
pause_time <- 3            # Number of seconds of waiting before going on with the script (this has an effect only if halt == TRUE)
n_points_for_arcsin <- 50  # Number of points on the arcsine pdf function
n_bins <- 20               # Number of bins in the histogram
rngseed <- 1111            # Seed for the random number generation

# User-intervention ends here

set.seed(rngseed)
coll_mat_simple <- matrix(data = NA, nrow = n_rep, ncol = 1)
coll_mat_combined <- matrix(data = NA, nrow = n_rep, ncol = n_toss-1)
colnames(coll_mat_simple) <- as.character(n_toss)
colnames(coll_mat_combined) <- c(2:n_toss-1)
for(k in 1:n_rep){
  rand_vec <- runif(n_toss, -1, 1)
  exp_vec <- character(0)
  step_vec <- integer(0)
  walk_vec <- 0
  avg_vec <- 0
  dist_vec <- integer(0)
  for(j in 1:length(rand_vec)){
    if(rand_vec[j] >= 0){
      #exp_vec <- c(exp_vec, "Heads")
      step_vec <- c(step_vec, -1)
      walk_vec <- c(walk_vec, walk_vec[length(walk_vec)] - 1)
      avg_vec <- c(avg_vec, walk_vec[length(walk_vec)] / j)
    } else {
      #exp_vec <- c(exp_vec, "Tails")
      step_vec <- c(step_vec, 1)
      walk_vec <- c(walk_vec, walk_vec[length(walk_vec)] + 1)
      avg_vec <- c(avg_vec, walk_vec[length(walk_vec)] / j)
    }
  }
  
  # Plotting and graphics unit starts here #
  
  if(intermed_plots != 0){
    if((n_rep %% intermed_plots) == 0){
      if(k%%intermed_plots == 0){
        print(k)
        walk_tab <- as.data.frame(cbind(c(0:n_toss), walk_vec, 10*avg_vec))
        colnames(walk_tab) <- c("toss","distance", "average")
        print(ggplot(data = walk_tab, aes(x = toss)) +
                geom_line(aes(y = distance, colour = "actual distance")) +
                geom_abline(slope = 0, intercept = 0, col = "red")+
                geom_line(aes(y = average, colour = "average"), size = 1.2) +
                labs(title = paste("Random walk representation of", as.character(n_toss), "coin tosses", sep = " "), color = "") +
                theme(plot.title = element_text(size = 12, hjust = 0.5), legend.position = "top") +
                scale_color_manual(values = c("average" = "green", "actual distance" = "black")))
        if (halt) {
          Sys.sleep(pause_time)
        }
      }
    }
  }
  coll_mat_simple[k,which(colnames(coll_mat_simple) == as.character(n_toss))]<- length(which(walk_vec > 0)) / length(walk_vec)
  for (i in 2:n_toss){
    coll_mat_combined[k,which(colnames(coll_mat_combined) == as.character(i))] <- length(which(walk_vec[2:i] > 0)) / length(walk_vec[2:i])
  }
}

for(i in 1:ncol(coll_mat_simple)){
  histogram_simple <- hist(coll_mat_simple[,i], plot = F, breaks = n_bins)
  x_hist_simple <- histogram_simple$mids
  dens_hist_simple <- histogram_simple$density
  hist_frame_simple <- as.data.frame(cbind(x_hist_simple, dens_hist_simple))
  arcsinpoints_simple <- matrix(data = 0, nrow = n_points_for_arcsin, ncol = 1)
  for (j in 1:n_points_for_arcsin){
    arcsinpoints_simple[j] <- 0 + (j/(n_points_for_arcsin+1))
  }
  darcsinpoints_simple <- darcsine(arcsinpoints_simple)
  data.frame(lapply(hist_frame_simple, "length<-", max(lengths(hist_frame_simple))))
  print(ggplot() + geom_col(aes(x = hist_frame_simple$x_hist_simple, y = hist_frame_simple$dens_hist_simple, fill = "Distance histogram"), size = 2, col = "green") +
          geom_line(aes(x = arcsinpoints_simple, y = darcsinpoints_simple, colour = "Arcsine pdf function"), size = 2) +
          labs(x = "The proportion of a run in which the distance is positive", y = "Frequency", title =paste("Histogram for", colnames(coll_mat_simple)[i], "coin tosses repeated", as.character(n_rep), "times", sep = " "), color = "", fill = "")+
          theme(plot.title = element_text(size = 12, hjust = 0.5), legend.position = "top")+
          scale_color_manual(values = c("Arcsine pdf function" = "blue")) + 
          scale_fill_manual(values = c("Histogram of the positive fraction of walks" = "green")))
}

histogram_combined <- hist(coll_mat_combined, plot = F, breaks = n_bins)
x_hist_combined <- histogram_combined$mids
dens_hist_combined <- histogram_combined$density
hist_frame_combined <- as.data.frame(cbind(x_hist_combined, dens_hist_combined))
arcsinpoints_combined <- matrix(data = 0, nrow = n_points_for_arcsin, ncol = 1)
for (j in 1:n_points_for_arcsin){
  arcsinpoints_combined[j] <- 0 + (j/(n_points_for_arcsin+1))
}
darcsinpoints_combined <- darcsine(arcsinpoints_combined)
data.frame(lapply(hist_frame_combined, "length<-", max(lengths(hist_frame_combined))))
print(ggplot() + geom_col(aes(x = hist_frame_combined$x_hist_combined, y = hist_frame_combined$dens_hist_combined, fill = "Distance histogram"), size = 2, col = "green") +
        geom_line(aes(x = arcsinpoints_combined, y = darcsinpoints_combined, colour = "Arcsine pdf function"), size = 2) +
        labs(x = "The proportion of a run in which the distance is positive", y = "Frequency", title =paste("Histogram for", as.character(n_toss), "coin tosses repeated", as.character(n_rep), "times", sep = " "), color = "", fill = "")+
        theme(plot.title = element_text(size = 12, hjust = 0.5), legend.position = "top")+
        scale_color_manual(values = c("Arcsine pdf function" = "blue")) + 
        scale_fill_manual(values = c("Histogram of the positive fraction of walks" = "green")))