metropolis1=function( n=1000 , rho=0.99 , start1=-1 , start2=-1 , tau=1 )
 {
  x1 <- c(1:n) 
  x2 <- c(1:n) 
  x1[1] <- start1
  x2[1] <- start2
  par( mfrow=c( 2,2 ) )
  for ( i in 2:n ) 
   {
    try <- x1[i-1] + rnorm( 1 , 0 , tau )
    temp <- exp( ( try^2 - x1[i-1]^2 - 2.0 * rho * x2[i-1] * ( try - x1[i-1] ) )/(-2.0*(1-rho^2)) ) 
    if ( runif(1) < temp ) 
     {
      x1[i] <- try
     } 
    else
     { 
      x1[i] <- x1[i-1]
     }
    try <- x2[i-1] + rnorm(1,0,tau )
    temp <- exp( ( try^2 - x2[i-1]^2 - 2.0 * rho * x1[i] *( try - x2[i-1]) )/(-2.0*(1-rho^2)) )
    if ( runif(1) < temp ) 
     {
      x2[i] <- try
     }
    else
     { 
      x2[i] <- x2[i-1]
     }
 

  plot( x1[1:i] , type="l" , xlab="Iteration" ,xlim=c(1,n),ylab="x1")
  plot( x2[1:i] , type="l" , xlab="Iteration",xlim=c(1,n),ylab="x2" )
  plot( x1[1:i] , x2[1:i] , type="p" , pch="." ,xlab="x1",ylab="x2")
  plot( x1[1:i] , x2[1:i] , type="l" ,col="red",xlab="x1",ylab="x2")
 }
 }