library(ts)

metunifnorm <- function(X0, n, scale)
{
  x <- X0
  output <- c(X0)
  for (i in 1:n)
  {
    y <- x + scale * (2*runif(1)-1)
    u <- runif(1)
    if ( (x^2 - y^2)/2 > log(u) )
      (x <- y)
    output <- c(output,x)
  }
  plot(output, type="l")
  output
}

metunifcauchy <- function(X0, n, scale)
{
  x <- X0
  output <- c(X0)
  for (i in 1:n)
  {
    y <- x + scale * (2*runif(1)-1)
    u <- runif(1)
    if ( (1+x^2)/(1+y^2) > u )
      (x <- y)
    output <- c(output,x)
  }
  plot(output, type="l")
  output
}

metnormnorm <- function(X0, n, scale)
{
  x <- X0
  output <- c(X0)
  for (i in 1:n)
  {
    y <- x + scale * rnorm(1)
    u <- runif(1)
    if ( (x^2 - y^2)/2 > log(u) )
      (x <- y)
    output <- c(output,x)
  }
  plot(output, type="l")
  output
}

metnormcauchy <- function(X0, n, scale)
{
  x <- X0
  output <- c(X0)
  for (i in 1:n)
  {
    y <- x + scale * rnorm(1)
    u <- runif(1)
    if ( (1+x^2)/(1+y^2) > u )
      (x <- y)
    output <- c(output,x)
  }
  plot(output, type="l")
  output
}


# e.g. x <- metgen( dnorm, rnorm, 0, 10000, 3)

metgen <- function(pi,proposal,X0, n, scale)
{
  output <- c(X0)
  x <- X0
  for (i in 1:n)
  {
    y <- x + scale * proposal(1)
    u <- runif(1)
    if ( log(pi(y)/pi(x)) > log(u) )
      (x <- y)
    output <- c(output,x)
  }
  plot(output, type="l")
  output
}

accept <- function(x)
{
  k<-length(x)
  number<-c(k-1)
  for (i in 1:(k-1))
  {
    if (x[i]==x[i+1])
      (number <- number-1)
  }
  number/(k-1)
}

indepsamp <- function(X0, n, lambda)
{
  output <- c(X0)
  x <- X0
  for (i in 1:n)
  {
    y <- rexp(1) / lambda 
    u <- runif(1)
    if ((x-y)*(1-lambda ) > log(u))
      (x<-y)
    output <- c(output,x)
  }
  plot(output, type="l")
  output
}

indepmeans <- function(m,lambda)
{
  themeans <- rep(0,m)
  for (i in 1:m)
    themeans[i] <- mean( indepsamp(1,10000,lambda) )
  themeans
  plot(density(themeans))
}

witch <- function(X0, n, w, delta)
{
  output <- c(X0)
  x <- X0
  for (i in 1:n)
  {
    y <- x + w * (2*runif(1)-1)
    if ( (y>0) && (y<1) && ( (y<delta) || (x>delta) ) )
      (x <- y)
    output <- c(output,x)
  }
  plot(output, type="l")
}