A (too) short introduction to R

In this first lab we will (if necessary) briefly see how one can use R.

0. Getting R (and a GUI)

First to use R, you need to install it! You can download it from here. Since the original GUI for R is a bit limited, you may want to install more advanced GUI. One popular option is RStudio. Once you have installed everything appropriately we can start using R!

1. Basic instructions

I just list here some simple computations.

x <- 1
y <- 3
x + y

[1] 4

x - y

[1] -2

x * y

[1] 3

x / y

[1] 0.3333333

x^y

[1] 1

## and many other stuffs: log, exp, sin, cos, gamma, ...

2. Vector, matrices, array, list

x <- c(1, 5, -3, 7, 8, -15)## vector of reals
x

[1]   1   5  -3   7   8 -15

y <- c("toto", "tutu", "tata")## vector of character strings
y

[1] "toto" "tutu" "tata"

A <- matrix(x, 3)##a matrix with 3 rows (hence 2 columns) which is "filled columnwise" with x
A

     [,1] [,2]
[1,]    1    7
[2,]    5    8
[3,]   -3  -15

B <- matrix(x, 3, byrow = TRUE)## filled rowwise
B

     [,1] [,2]
[1,]    1    5
[2,]   -3    7
[3,]    8  -15

C <- matrix(y, 3, 4)## recycle the vector y to fill in the matrix
C

     [,1]   [,2]   [,3]   [,4]  
[1,] "toto" "toto" "toto" "toto"
[2,] "tutu" "tutu" "tutu" "tutu"
[3,] "tata" "tata" "tata" "tata"

## Tensors are also available from the array function
D <- array(1:60, dim = c(3, 4, 5))## A 3d tensor
E <- array(1:60, dim = c(5, 2, 3, 2))## A 4d tensor
D

, , 1

     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12

, , 2

     [,1] [,2] [,3] [,4]
[1,]   13   16   19   22
[2,]   14   17   20   23
[3,]   15   18   21   24

, , 3

     [,1] [,2] [,3] [,4]
[1,]   25   28   31   34
[2,]   26   29   32   35
[3,]   27   30   33   36

, , 4

     [,1] [,2] [,3] [,4]
[1,]   37   40   43   46
[2,]   38   41   44   47
[3,]   39   42   45   48

, , 5

     [,1] [,2] [,3] [,4]
[1,]   49   52   55   58
[2,]   50   53   56   59
[3,]   51   54   57   60

l <- list(component1 = x, component2 = C)## list the most "powerfull type" can mix different object class within a list
l

$component1
[1]   1   5  -3   7   8 -15

$component2
     [,1]   [,2]   [,3]   [,4]  
[1,] "toto" "toto" "toto" "toto"
[2,] "tutu" "tutu" "tutu" "tutu"
[3,] "tata" "tata" "tata" "tata"

D <- cbind(x, y)## column binding to create a 2 column matrix (of which class???)
D

     x     y     
[1,] "1"   "toto"
[2,] "5"   "tutu"
[3,] "-3"  "tata"
[4,] "7"   "toto"
[5,] "8"   "tutu"
[6,] "-15" "tata"

D2 <- rbind(x, y)## row binding
D2

  [,1]   [,2]   [,3]   [,4]   [,5]   [,6]  
x "1"    "5"    "-3"   "7"    "8"    "-15" 
y "toto" "tutu" "tata" "toto" "tutu" "tata"

3. Linear algebra

x <- 1:3##integers from 1 to 3 (included)
y <- seq(2.1, 15, length = 3)##evenly space points between 2.1 ane 15 (total length 3)

x + y

[1]  3.10 10.55 18.00

x * y ## componentwise multiplication

[1]  2.1 17.1 45.0

x %*% y ## inner product

     [,1]
[1,] 64.2

A <- diag(5, 3)## 3 x 3 diagonal matrix with 5 on the diagonal
B <- matrix(1:9, 3)
t(B) ## Matrix transpose

     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6
[3,]    7    8    9

A %*%  ## matrix multiplication
A * A ## componentwise multiplication (Hadamard product)

     [,1] [,2] [,3]
[1,]  125    0    0
[2,]    0  125    0
[3,]    0    0  125

solve(A)## matrix inverse

     [,1] [,2] [,3]
[1,]  0.2  0.0  0.0
[2,]  0.0  0.2  0.0
[3,]  0.0  0.0  0.2

eigen(A)## eigen decomposition

eigen() decomposition
$values
[1] 5 5 5

$vectors
     [,1] [,2] [,3]
[1,]    0    0    1
[2,]    0    1    0
[3,]    1    0    0

svd(A)## singular value decomposition

$d
[1] 5 5 5

$u
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1

$v
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1

4. Random number generation

rnorm(5, 0, 4)##5 independent draws from a gaussian distribution with mean 0 and *standard deviation* 4.

[1]  2.0176446  1.3996484 -2.3870872  0.1009725  3.8961008

runif(5, -4, 3)## 5 independent draws from a U(-4, 3) distribution

[1] -1.5098755  2.2585768 -1.0820546 -2.7980153  0.8389459

5. Basic plots

## Generate a random walk
randomwalk <- cumsum(rnorm(1000))
plot(randomwalk)

plot(randomwalk, type = "l")

plot(randomwalk, type = "o")

x <- seq(0, pi, length = 5)
y <- sin(x)
plot(x, y, xlab = "x", ylab = "sin(x)")

plot(sin, from = 0, to = 2 * pi, col = "orange", lwd = 2)##lwd for Line WiDth

x <- rnorm(100)
hist(x, freq = FALSE)##freq = FALSE so that area = 1
plot(dnorm, from = -4, to = 4, col = "green", lwd = 2, lty = 2, add = TRUE)## add = TRUE to add the plot to the current graphic

x <- rnorm(100)
y <- rnorm(100, 1)
boxplot(cbind(x, y))##boxplots

6. Importing / Exporting data

Well read the documentation of the read.table, write.table and save function.

7. For, if, else, loops

for (i in 1:3)
  print(i)

[1] 1
[1] 2
[1] 3

astring <- "what"

if (astring == "who"){
  outstring <- "what is"
  a <- 1
} else if (astring == "which"){
  outstring <- "what?"
  a <- 2
} else {
  outstring <- "is what?"
  a <- 3
}

print(outstring)

[1] "is what?"

8. Creating function

myfunction <- function(arg1, arg2, optional.argument = "take log"){
  output <- arg1 + arg2
  
  if (optional.argument == "take log")
    output <- log(output)
  
  return(output)
}

myfunction(2, 3)

[1] 1.609438

myfunction(2, 3, "do not take log")

[1] 5

Watch out an R function can only return a single object so if you want to ouput several quantities you will have to store them within a list, i.e., as below

myfunction2 <- function(arg1, arg2, optional.argument = "take log"){
  output <- arg1 + arg2
  
  if (optional.argument == "take log")
    output <- log(output)
  
  return(list(output = output, message = optional.argument))
}

myfunction2(2, 3)

$output
[1] 1.609438

$message
[1] "take log"

9. Numerical optimization

rosenbrock <- function(x) {   ## Rosenbrock Banana function
    x1 <- x[1]
    x2 <- x[2]
    100 * (x2 - x1 * x1)^2 + (1 - x1)^2
}

n.grid <- 100
x <- seq(-2, 2, length = n.grid)
y <- seq(-3, 3, length = n.grid)
grid <- expand.grid(x, y)
z <- apply(grid, 1, rosenbrock)##compute the function on the grid
z <- matrix(z, n.grid)

## plot of it (just take care of the persp line of code not the previous ones)
jet.colors <- colorRampPalette( c("blue", "green") )
nbcol <- 250
color <- jet.colors(nbcol)
# Compute the z-value at the facet centres
zfacet <- z[-1, -1] + z[-1, -n.grid] + z[-n.grid, -1] + z[-n.grid, -n.grid]
# Recode facet z-values into color indices
facetcol <- cut(zfacet, nbcol)
persp(x, y, z, col = color[facetcol], phi = 40, theta = -30, border = NA)

## minimize the banana function
init <- c(0, 0)## minimum is at (1, 1)
optim(init, rosenbrock)

$par
[1] 0.9999564 0.9999085

$value
[1] 3.729052e-09

$counts
function gradient 
     169       NA 

$convergence
[1] 0

$message
NULL

optim(init, rosenbrock, method = "BFGS")## changing optimizer

$par
[1] 0.9998000 0.9996001

$value
[1] 3.998081e-08

$counts
function gradient 
      63       26 

$convergence
[1] 0

$message
NULL

nlm(rosenbrock, init)## another (the one I used often) optimizer

$minimum
[1] 4.023727e-12

$estimate
[1] 0.999998 0.999996

$gradient
[1] -7.328280e-07  3.605689e-07

$code
[1] 1

$iterations
[1] 26