Untitled

Import some vis packages

library(ggplot2)
library(tidyverse)

## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --

## v tibble  3.1.0     v dplyr   1.0.5
## v tidyr   1.1.3     v stringr 1.4.0
## v readr   1.4.0     v forcats 0.5.1
## v purrr   0.3.4

## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

# Given function Q4.1. Call it f
f <- function(x){
  
  return(sin(x))
}
g <- function(x){
  
  return(cos(x))
}
g_prime <- function(x){
  
  return(-sin(x))
}

# Since we use limits, we need a package called "Ryacas" - install if you don't already have it
library(Ryacas)     # version 1.1.1

## 
## Attaching package: 'Ryacas'

## The following object is masked from 'package:purrr':
## 
##     simplify

## The following object is masked from 'package:stats':
## 
##     integrate

## The following objects are masked from 'package:base':
## 
##     %*%, diag, diag<-, lower.tri, upper.tri

x <- ysym("x")
a <- ysym("a")
b <- ysym("b")
h <- ysym("h")
y <- ysym("y")
# define x as a symbolic variable
#lim(sin(x)/x, x, 0) # limit of sin(x)/x as x approaches 0
## 1

first_derivative <- function(funct,a,h){
  return(lim((funct(a+h) - funct(a))/h, h, 0))
}

second_derivative <- function(funct,a,h){
  return(lim((funct(a+h) - 2*funct(a)+f(a-h))/h^2, h, 0))
}
print("First derivative:")

## [1] "First derivative:"

first_derivative(f, x, h)

## y: Cos(x)

# Hard code the first derivative:
g <- function(x){
  
  return(cos(x))
}
print("Evaluated at x=5:")

## [1] "Evaluated at x=5:"

g(5)

## [1] 0.2836622

print("#######################################")

## [1] "#######################################"

print("Second derivative:")

## [1] "Second derivative:"

first_derivative(g, x, h) # this is just a test

## y: -Sin(x)

second_derivative(f, x, h) 

## y: (-2*Sin(x))/2

print("Evaluated at x=5:")

## [1] "Evaluated at x=5:"

g_prime(5)

## [1] 0.9589243

# Visualise
ggplot(data = data.frame(x = 0), mapping = aes(x = x)) + 
  stat_function(fun = f, size = 1.25, alpha = 0.75, colour= "blue")  +
  xlim(-6,6)

ggplot(data = data.frame(x = 0), mapping = aes(x = x)) + 
  stat_function(fun = g, size = 1.25, alpha = 0.75, colour= "red")  +
  xlim(-6,6)+ggtitle("First derivative of f(x)=sin(x) is cos(x)")

ggplot(data = data.frame(x = 0), mapping = aes(x = x)) + 
  stat_function(fun = g_prime, size = 1.25, alpha = 0.75, colour= "green")  +
  xlim(-6,6) + ggtitle("Second derivative of f(x)=sin(x) is -sin(x)")

4.2

# The full function:
z_func <- function(x,y){
  return((3*x*y) + (x^2) + (y^2))
}

# When x = 5
#zy <- function(y){
# return(15*y + 25 + y^2)
#}

# when y = 3
#zx <- function(x){
#  return(9*x + 9 + x^2)
#}

# Graph the function in full to see:
##################################
x <- seq(-10, 10, length= 30)
y <- x
z <- outer(x, y, z_func)
z[is.na(z)] <- 1
op <- par(bg = "white")
persp(x, y, z, theta = 30, phi = 30, expand = 0.5, col = "lightblue")

###################################


# This will give a more colourful view of the same graph:
x <- seq(-10, 10, length= 30)
y <- x
z <- outer(x, y, z_func)
z[is.na(z)] <- 1

require(lattice)

## Loading required package: lattice

wireframe(z, drape=T, col.regions=rainbow(100))

# Now that we know what it looks like, we can determine the derivatives:
# Since the function is in the 3D-plane, its derivatives should be in the 2D plane.

# Let's check the answer by using the built-in function D:
f<-expression(x^2+y^2+3*x*y)
first_wrt_x <- D(f,'x')
first_wrt_x  #Q4.2.1

## 2 * x + 3 * y

 # at x=5 and y=3, this is equal to dz/dx = 19
first_wrt_y <- D(f,'y')
first_wrt_y #Q4.2.2

## 2 * y + 3 * x

 # at x=5 and y=3, this is equal to dz/dy = 21

# Now for the second derivatives:
second_wrt_x_and_y <- D(D(f, 'x'), 'y')
second_wrt_x_and_y #Q4.2.3

## [1] 3

 # at x=5 and y=3, this is equal to dz^2/dxdy = 3

second_wrt_x <- D(D(f, 'x'), 'x')
second_wrt_x #Q4.2.4

## [1] 2

 # at x=5 and y=3, this is equal to dz^2/dx^2 = 2

# Since we use limits, we need a package called "Ryacas" - install if you don't already have it
library(Ryacas)     # version 1.1.1
x <- ysym("x")      # define x as a symbolic variable
# lim(sin(x)/x, x, 0) # limit of sin(x)/x as x approaches 0
## 1

a <- ysym("a")
b <- ysym("b")
h <- ysym("h")
y <- ysym("y")

partial_derivative_wrt_x <- function(funct,a,b,h){
  return(lim((funct(a+h,b) - funct(a,b))/h, h, 0))
}

partial_derivative_wrt_y <- function(funct,a,b,h){
  return(lim((funct(a,b+h) - funct(a,b))/h, h, 0))
}
z_func <- function(x,y){
  return((3*x*y) + (x^2) + (y^2))
}

# Q4.2.1
print("Question 4.2 - first bullet")

## [1] "Question 4.2 - first bullet"

partial_derivative_wrt_x(z_func, x, y, h)

## y: 3*y+2*x

first_deriv_wrt_x <- function(x,y){
  return((3*y) + (2*x))
}
# Hard code for x=5 and y=3:
print("Evaluated at x=5 and y=3:")

## [1] "Evaluated at x=5 and y=3:"

first_deriv_wrt_x(5,3)

## [1] 19

print("#######################################")

## [1] "#######################################"

# Q4.2.2:
print("Question 4.2 - second bullet")

## [1] "Question 4.2 - second bullet"

partial_derivative_wrt_y(z_func, x, y, h)

## y: 3*x+2*y

first_deriv_wrt_y <- function(x,y){
  return((3*x) + (2*y))
}
# Hard code for x=5 and y=3:
print("Evaluated at x=5 and y=3:")

## [1] "Evaluated at x=5 and y=3:"

first_deriv_wrt_y(5,3)

## [1] 21

print("#######################################")

## [1] "#######################################"

# Q4.2.3:
print("Question 4.2 - third bullet")

## [1] "Question 4.2 - third bullet"

partial_derivative_wrt_y(first_deriv_wrt_x, x, y, h)

## y: 3

print("#######################################")

## [1] "#######################################"

# Q4.2.4:
print("Question 4.2 - forth bullet")

## [1] "Question 4.2 - forth bullet"

partial_derivative_wrt_x(first_deriv_wrt_x, x, y, h)

## y: 2