%\documentclass[handout]{beamer} \documentclass{beamer} \usetheme{Rochester} \usecolortheme{wolverine} %\definecolor{UBCblue}{rgb}{0.04706, 0.13725, 0.26667} % UBC Blue (primary) %\usecolortheme[named=UBCblue]{structure} \definecolor{LEEgreen}{rgb}{0.04,0.667,0.137} \usepackage{amsmath} \usepackage{multirow} \usepackage{latexsym} \usepackage{amssymb} \usepackage{pifont} %\usepackage{epic} \usepackage{graphicx} \usepackage{skull} %\usepackage{multicol} %\usepackage{cancel} %\usepackage{pdfpages} \usepackage{listings} \usepackage{mathtools} \usepackage{fancybox} \usepackage{ulem} \DeclarePairedDelimiter\abs{\lvert}{\rvert} \addtolength{\parskip}{1ex} \setbeamercovered{transparent=5} \definecolor{ForestGreen}{cmyk}{0.5, 0, 0.5, 0.5} \definecolor{DarkBlue}{rgb}{0.0, 0.08, 0.45} \definecolor{BrickRed}{cmyk}{0, .89, .94, .28} \definecolor{UBCblue}{rgb}{0.04706, 0.13725, 0.66667} % UBC Blue (primary) \def \hi{\textcolor{BrickRed}} \def \E{{\textrm{E}}} \newcommand{\ex}[2][]{\ensuremath{\mathrm{E}_{#1}\!\left[#2\right]}} \newcommand{\var}[2][]{\ensuremath{\mathrm{Var}_{#1}\left(#2\right)}} \newcommand{\cov}[2][]{\ensuremath{\mathrm{Cov}\left(#2\right)}} \newcommand{\cor}[2][]{\ensuremath{\mathrm{Cor}\left(#2\right)}} \newcommand{\given}{\ensuremath{\;|\;}} \def \R{\mathbb{R}} \def \proof{\noindent{\bf Proof: }} \def \eproof{\qed\vskip .3cm} \def \ss{S} \def \union{\cup} \def \intersect{\cap} \newcommand{\vsp}{\vspace{0.3cm}} \newcommand{\vspl}{\vspace{0.15cm}} \newenvironment{answers}{\color{BrickRed} ~}{\color{black} ~} \newcommand{\pr}[2][]{\ensuremath{\mathrm{Pr}_{#1}\!\left(#2\right)}} \newcommand{\perm}[2]{\operatorname{P}^{#1}_{#2}} \newcommand{\sol}{\setbeamercolor{structure}{fg=green!60!black}} \renewcommand{\vec}[1]{\underline{#1}} \title{Data Analytics 1: Introduction to R and RStudio} \author[Lee Fawcett \& Joe Matthews]{Dr. Lee Fawcett \& Dr. Joe Matthews} \date{\it NUiS Overseas Immersion Programme 2026} %\usepackage[usenames, dvipsnames]{xcolor} \begin{document} %%From web for code/highlighting \definecolor{mygreen}{rgb}{0,0.6,0} \definecolor{mygray}{rgb}{0.5,0.5,0.5} \definecolor{mymauve}{rgb}{0.58,0,0.82} \definecolor{light-gray}{gray}{0.95} \lstset{ backgroundcolor=\color{light-gray}, % choose the background color; you must add \usepackage{color} or \usepackage{xcolor}; should come as last argument basicstyle=\footnotesize, % the size of the fonts that are used for the code breakatwhitespace=false, % sets if automatic breaks should only happen at whitespace breaklines=true, % sets automatic line breaking captionpos=b, % sets the caption-position to bottom commentstyle=\color{mygreen}, % comment style deletekeywords={...}, % if you want to delete keywords from the given language escapeinside={\%*}{*)}, % if you want to add LaTeX within your code extendedchars=false, % lets you use non-ASCII characters; for 8-bits encodings only, does not work with UTF-8 frame=single, % adds a frame around the code keepspaces=true, % keeps spaces in text, useful for keeping indentation of code (possibly needs columns=flexible) keywordstyle=\color{blue}, % keyword style language=R, % the language of the code morekeywords={*,...}, % if you want to add more keywords to the set numbers=left, % where to put the line-numbers; possible values are (none, left, right) numbersep=5pt, % how far the line-numbers are from the code numberstyle=\tiny\color{mygray}, % the style that is used for the line-numbers rulecolor=\color{black}, % if not set, the frame-color may be changed on line-breaks within not-black text (e.g. comments (green here)) showspaces=false, % show spaces everywhere adding particular underscores; it overrides 'showstringspaces' showstringspaces=false, % underline spaces within strings only showtabs=false, % show tabs within strings adding particular underscores stepnumber=1, % the step between two line-numbers. If it's 1, each line will be numbered stringstyle=\color{mymauve}, % string literal style tabsize=2, % sets default tabsize to 2 spaces title=\lstname % show the filename of files included with \lstinputlisting; also try caption instead of title } \begin{frame} \frametitle{\small \textit{\color{purple}{2pm, Wednesday 8th July}}} \titlepage \end{frame} \frame{ \frametitle{Who we are} Lee Fawcett \begin{itemize} \item Reader in Statistics and Associate Dean (Global) \pause \item Teach: Computational Statistics, Bayesian Statistics, Extreme Value Theory, `Global' \pause \item Research: Environmental Extremes, Statistics and Road Safety, Stats Education \pause \end{itemize} \vsp Joe Matthews \begin{itemize} \item Senior Lecturer in Statistical Data Science \pause \item Teach: Data Science (MSc programmes), `Global' \pause \item Research: Statistics and Road Safety, Stats Education \end{itemize} } \frame{ \frametitle{This short course} \begin{enumerate} \item Introdution to R \& RStudio \pause \item Data Visualisation\pause \item Multivariate Data Analysis \pause \item We'll see...! \pause \end{enumerate} In all sessions: A mix of interactive lectures and practical sessions! } \subsection{Books and other resources} \frame{ \frametitle{Books and other resources} Some useful books include: \pause \begin{itemize} \item W.N. Venables and D.M. Smith: {\it An Introduction To R.} Free to download at: {\small \begin{center} {\color{blue} \url{https://cran.r-project.org/doc/manuals/R-intro.pdf}} \end{center}} \item H. Wickham, M Çetinkaya-Rundel and G. Grolemund: {\it R for Data Science.} \item A. Field: {\it Discovering Statistics Using R \& RStudio.} \end{itemize}} \subsection{General background} \frame{ \frametitle{General background} In this class we will learn about computational tools for performing statistical analysis. \vsp We will use the R programming language. \vsp R is very widely used in universities, businesses and other organizations due to its many advantages over other similar software packages. \begin{itemize} \item It is available for \textcolor{blue}{\bf free} \item It is `\textcolor{blue}{\bf open source}' and can be used for a very wide range of statistical and scientific applications \item \textcolor{blue}{\bf Libraries} of R functions and data-structures have been created and are readily available on the web \end{itemize} } \frame{ \frametitle{General background} Learning to use R is important because: \begin{itemize} \item [(i)] It is the foremost tool in \textcolor{blue}{\bf modern computational statistics} \item [(ii)] It is a great way to learn general concepts in \textcolor{blue}{\bf programming} which can be used with other languages or in other contexts \item [(iii)] It can be used to illustrate ideas in \textcolor{blue}{\bf Statistics \& Data Science} \end{itemize} } \frame{ \frametitle{Preliminaries} R is installed on all University machines. \pause \vsp In this module we will be using \textcolor{blue}{\bf RStudio}, which is an R \textit{Integrated Development Environment}, or IDE. \vsp This just makes R much easier to use and navigate, through the four windows that you see when you open RStudio (text editor, console window, file explorer and graphics display). } \frame{ \frametitle{Preliminaries} In today's brief intro, we will quickly review: \begin{itemize} \item [(i)] Data types and vectors in R \item [(ii)] Dataframe manipulation \item [(iii)] Data summaries \item [(iv)] Functions and control statements \end{itemize} \pause But first... } \begin{frame}[fragile] \frametitle{Using R like a calculator} At its most basic level, R can be used like a big calculator: \pause \begin{lstlisting} #%* \textcolor{ForestGreen}{Let's do addition:} *) > 174+895 [1] 1069 \end{lstlisting} \pause \vspace{-0.5cm} \begin{lstlisting} #%* \textcolor{ForestGreen}{Now try division:} *) > 96/12 [1] 8 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Using R like a calculator} \begin{lstlisting} #%* \textcolor{ForestGreen}{Notice the command prompt $>$} *) #%* \textcolor{ForestGreen}{(I'll not always include that in the slides)} *) #%* \textcolor{ForestGreen}{Also notice that nothing following \# is evaluated} *) \end{lstlisting} \pause \vspace{-0.5cm} \begin{lstlisting} #%* \textcolor{ForestGreen}{Trig functions? No problem:} *) > sin(pi/2) [1] 1 \end{lstlisting} \pause \vspace{-0.5cm} \begin{lstlisting} #%* \textcolor{ForestGreen}{Powers? Sorted.} *) > sqrt(12^2) [1] 12 \end{lstlisting} \pause \vspace{-0.75cm} R can do basic (and not-so-basic) computations on vectors and matrices too! \end{frame} \frame{ \frametitle{The Movies dataset} I like watching movies -- most people do -- so I often use the IMDB movies database to demonstrate key concepts in RStudio (don't worry, there will be more relevant examples for you soon!). \vsp IMDB is a website devoted to collecting movie data supplied by studios and fans. \vsp More information can be found online: {\color{blue} \begin{center} \url{http://imdb.com/help/show_leaf?about} \end{center} } } \frame{ \frametitle{The Movies dataset} \input movietableforslides } \begin{frame}[fragile] \frametitle{The Movies dataset} I have written an add-on R package which will allows easy access to many interesting datasets (including the movies dataset, but others you might be interested in too!). Installing the package is straightforward: % Turn off + signs in R listings % Other sweave tricks %\DefineVerbatimEnvironment{Sinput}{Verbatim} {xleftmargin=2em} %\DefineVerbatimEnvironment{Soutput}{Verbatim}{xleftmargin=2em} %\DefineVerbatimEnvironment{Scode}{Verbatim}{xleftmargin=2em} %\fvset{listparameters={\setlength{\topsep}{0pt}}} %\renewenvironment{Schunk}{\vspace{\topsep}}{\vspace{\topsep}} \begin{lstlisting} install.packages("ncldata", repos = "http://R-Forge.R-project.org", type = "source") \end{lstlisting} \pause \vspace{-0.3cm} To load the package, type: \begin{lstlisting} library(ncldata) \end{lstlisting} \vspace{-0.3cm} Then, to specifically load the movies dataset, type: \begin{lstlisting} data(movies) \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{The Movies dataset} Alternatively, you can by-pass the package altogether as I've posted the dataset to my personal webpage: \begin{lstlisting} movies = read.table("https://www.mas.ncl.ac.uk/~nlf8/basicR/movies.txt", header=TRUE) attach(movies) \end{lstlisting} \pause \vspace{-0.5cm} OK, now on to the meat of today's session...! \end{frame} \subsection{Data types in R} \begin{frame}[fragile] \frametitle{Data types in R} R has a variety of data types. For example: \begin{lstlisting} v = TRUE #%* \textcolor{ForestGreen}{`Logical' or `boolean' - TRUE or FALSE} *) w = "fred" #%* \textcolor{ForestGreen}{`Character' or `string'} *) x = 5.0 #%* \textcolor{ForestGreen}{`Double' or `numerical'} *) \end{lstlisting} \vspace{-0.5cm} There are also some `special' data types: \begin{lstlisting} y = 5/0 #%* \textcolor{ForestGreen}{Returns \texttt{Inf} }*) z = y-y #%* \textcolor{ForestGreen}{Returns \texttt{NaN} }*) \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Data types in R} Another important data type is \texttt{NA}, short for `Not Available' and used for missing values. \vsp A summary of the data types in R is given in the Table below. \vsp {\small %\begin{table}[!h] \begin{tabular}{l|l l l l} \hline Type &Example 1 &Example 2 &Example 3 &Example 4 \\ \hline Doubles &2 &3.1242 &-45.6 &4e-10 \\ Logicals/Boolean&TRUE &FALSE & & \\ Characters &``Fred'' &``TRUE'' &``Male'' &``2'' \\ Infinity &Inf &5/0 & & \\ Other &NaN &NA & & \\ \hline \end{tabular} %\caption{A summary of some data types in R.} %\label{table:datatypes} %\end{table} } \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#1}} In the code below, what is the value of x? \begin{lstlisting} x = 5 x = 2*x + 1 x + 5 \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 5 \item [2.] 11 \item [3.] 16 \item [4.] 22 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#1}} In the code below, what is the value of x? \begin{lstlisting} x = 5 x = 2*x + 1 x + 5 \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 5 \item [2.] \textcolor{red}{\bf 11} \item [3.] 16 \item [4.] 22 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#2}} In the code below, what type of variable is x? \begin{lstlisting} x = 5 x = 2*x + 1 x + 5 \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] Character \item [2.] Double \item [3.] Logical \item [4.] NA \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#2}} In the code below, what type of variable is x? \begin{lstlisting} x = 5 x = 2*x + 1 x + 5 \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] Character \item [2.] \textcolor{red}{\bf Double} \item [3.] Logical \item [4.] NA \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#3}} In the code below, what type of variable is y? \begin{lstlisting} x = 0 y = 2/0 y = "Fred" \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] NaN \item [2.] Double \item [3.] NA \item [4.] Character \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#3}} In the code below, what type of variable is y? \begin{lstlisting} x = 0 y = 2/0 y = "Fred" \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] NaN \item [2.] Double \item [3.] NA \item [4.] \textcolor{red}{\bf Character} \end{itemize} \end{frame} \subsection{Vectors} \begin{frame}[fragile] \frametitle{Vectors} Vectors are the most basic of all data structures, but are used in almost all R code. \vsp An R vector contains $n$ values of the same type, where $n$ can be zero. \vsp For example: \begin{lstlisting} c(0, 1, 2, 3, 4, 5) [1] 0 1 2 3 4 5 (first_vec = c(0, 1, 2, 3, 4, 5)) [1] 0 1 2 3 4 5 (second_vec = c("Male", "Female", "Male")) [1] "Male" "Female" "Male" x = 5.2 y = "Fred" \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Vectors} There are special functions in R that can be used to determine the variable type. For example: \begin{lstlisting} is.double(x) [1] TRUE is.character(x) [1] FALSE is.character(y) [1] TRUE is.vector(first_vec) [1] TRUE \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Vectors} And of course, we can use the \texttt{length} command to determine the length of a vector: \begin{lstlisting} length(first_vec) [1] 6 length(second_vec) [1] 3 \end{lstlisting} \pause And you can use the \texttt{seq} command to generate sequence vectors: \begin{lstlisting} (x1 = seq(1, 6)) [1] 1 2 3 4 5 6 (x2 = seq(-4, 4, by=2)) [1] -4 -2 0 2 4 (x3 = seq(3, -2, by=-0.5)) [1] 3.0 2.5 2.0 1.5 1.0 0.5 0.0 -0.5 -1.0 -1.5 [11] -2.0 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Vectors} Other useful commands for vectors: \begin{lstlisting} rev(x1) #%* \textcolor{ForestGreen}{Reverses the order of the vector}*) [1] 6 5 4 3 2 1 sort(x3) #%* \textcolor{ForestGreen}{Sorts the vector into ascending order}*) [1] -2.0 -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 sort(second_vec) #%* \textcolor{ForestGreen}{Can be applied to characters: alphabetical order}*) [1] "Female" "Male" "Male" sum(x2); sum(x3); sum(second_vec) [1] 0 [1] 5.5 Error in sum(second_vec) : invalid 'type' (character) of argument unique(second_vec) #%* \textcolor{ForestGreen}{Lists all unique elements of the vector}*) [1] "Male" "Female" which(x2==2) #%* \textcolor{ForestGreen}{Gives index/indices of particular elements}*) [1] 4 \end{lstlisting} \end{frame} \subsubsection{Vector operations} \begin{frame}[fragile] \frametitle{Vector operations} When our data are in a vector structure, we can apply standard operations to the entire vector. For example: \begin{lstlisting} (x = seq(-4, 4)) [1] -4 -3 -2 -1 0 1 2 3 4 x*x [1] 16 9 4 2 1 0 1 2 4 9 16 x-5 [1] -9 -8 -7 -6 -5 -4 -3 -2 -1 x + x [1] -8 -6 -4 -2 0 2 4 6 8 \end{lstlisting} \end{frame} \subsubsection{Extracting elements from vectors}\label{section:vectorextraction} \begin{frame}[fragile] \frametitle{Extracting elements from vectors} R has a number of useful methods that we can use to extract subsets of our data. \vsp For example, square brackets: \begin{lstlisting} first_vec[2] [1] 1 second_vec[2:3] [1] "Female" "Male" first_vec[4:2] [1] 3 2 1 (z = 1:10) [1] 1 2 3 4 5 6 7 8 9 10 z[z<3 | z>=7] #%* \textcolor{ForestGreen}{Values of z that are less than 3 \textbf{or} at least 7}*) [1] 1 2 7 8 9 10 \end{lstlisting} \vspace{-0.5cm} \pause We will look at extraction from two-dimensional R objects shortly. \end{frame} \begin{frame}[fragile] \frametitle{Extracting elements from vectors} Another useful way of extracting elements from vectors is to use logical vectors. \vsp Recall that the logical data type in R can take only the values \texttt{TRUE} and \texttt{FALSE}. Hence: \begin{lstlisting} A = TRUE B = FALSE !A #%* \textcolor{ForestGreen}{Read as ``Not A"}*) [1] FALSE !B [1] TRUE A & B #%* \textcolor{ForestGreen}{Need ``both A \textbf{and} B"}*) [1] FALSE A | B #%* \textcolor{ForestGreen}{Need ``either A \textbf{or} B"}*) [1] TRUE \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Extracting elements from vectors} Then we could declare a vector of such logicals to help us extract from another vector: \begin{lstlisting} (logic1 = c(TRUE, FALSE, TRUE, FALSE)) [1] TRUE FALSE TRUE FALSE (vec1 = seq(3, 6)) [1] 3 4 5 6 vec1[logic1] #%* \textcolor{ForestGreen}{Will extract the elements in vec1 where logic1 is true}*) [1] 3 5 \end{lstlisting} \pause We can also sum a logical vector. For example: \begin{lstlisting} sum(logic1) [1] 2 sum(vec1[logic1]) [1] 8 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Example (One for later...)} \includegraphics[width=1cm]{ironman.jpg} In RStudio, commit the following lines of code. Both lines read in files stored on a webpage. \textbf{Do not copy-and-paste this code -- you must type it in carefully!} \begin{lstlisting} x = scan("http://www.mas.ncl.ac.uk/~nlf8/basicR/x.txt") y = as.logical(scan("http://www.mas.ncl.ac.uk/~nlf8/basicR/y.txt", what="logical")) \end{lstlisting} \vspace{-1cm} {\tiny \begin{center} \begin{tabular}{l l} Questions about x &Questions about y \\ \hline (a) What type of data is in the vector x? &(a) What type of data is in the vector y? \\ (b) How many values are in the vector x? &(b) What is the first value of y?\\ (c) What is the smallest value in the vector x? &(c) What is the last value of y?\\ (d) What is the value of x in position 3500? &(d) How many times does TRUE occur? \\ (e) What is the 97th smallest value in x? &(e) When is the 472nd occurrence of TRUE? \\ (f) Find $\sum_{i=5612}^{8664}\text{x}_{i}$. &(f) How many occurrences of TRUE are there \\ & \textcolor{white}{LL} between positions 504 and 593 (inclusive)? \\ \hline \end{tabular} \end{center} } \end{frame} \subsection{Dataframes}\label{sec:dataframes} \begin{frame} \frametitle{Dataframes} A dataframe is a special kind of object in R. \vsp We use dataframes for storing and managing datasets that have a rectangular structure. \vsp Typically the rows correspond to \textcolor{blue}{\bf cases} and the columns \textcolor{blue}{\bf variables}. \vsp The crucial difference between a dataframe and a matrix is that all values in a matrix must be of the same type, whereas in a dataframe we can have responses for cases on variables of different types. \end{frame} \begin{frame} \frametitle{Dataframes} For example, we might have survey results on 10 individuals, who have each been asked 5 questions. \begin{itemize} \item This would give a dataframe with 10 rows and 5 columns \item We call this a ``\textcolor{blue}{\bf 10 by 5 dataframe}" (generally, we refer to $r \times c$ dataframes) \item The variable types need not all be the same: \begin{itemize} \item [--] The first column could be a vector of doubles recording the individuals' height in metres \item [--] The second column could be a vector of logicals taking the value \texttt{TRUE} if they like Marmite \item [--] The third column could be a vector of characters taking the name of their favourite book, etc. \end{itemize} \end{itemize} \end{frame} \subsubsection{Creating and querying your own dataframe} \begin{frame}[fragile] \frametitle{Creating and querying your own dataframe} You can construct your own dataframe using R's \texttt{data.frame} command. For example: \begin{lstlisting} height_m = c(1.81, 1.75, 1.58) weight_kg = c(NA, NA, 70.3) marmite = c(TRUE, FALSE, FALSE) book = c("Midnight in Chernobyl", "To Kill a Mockingbird", "Life of Pi") (small.survey = data.frame(height_m, weight_kg, marmite, book)) height_m weight_kg marmite book 1 1.81 NA TRUE Midnight in Chernobyl 2 1.75 NA FALSE To Kill a Mockingbird 3 1.58 70.3 FALSE Life of Pi \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Creating and querying your own dataframe} If you had a large dataframe and wanted to quickly know its size, you could use the \texttt{dim} command. For example: \begin{lstlisting} dim(small.survey) [1] 3 4 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Creating and querying your own dataframe} We can extract elements from a dataframe using square brackets (just as we demonstrated earlier for vector extraction). \vsp If the columns have labels, we can also use the \$ command. \begin{lstlisting} small.survey[2, 4] #%* \textcolor{ForestGreen}{Extract row 2, column 4}*) [1] To Kill a Mockingbird small.survey[ , 4] #%* \textcolor{ForestGreen}{Extract all rows from column 4}*) [1] Midnight in Chernobyl To Kill a Mockingbird Life of Pi small.survey$book #%* \textcolor{ForestGreen}{Does exactly the same as line 4!}*) [1] Midnight in Chernobyl To Kill a Mockingbird Life of Pi small.survey[3, ] #%* \textcolor{ForestGreen}{Extract all columns from row 3}*) height_m weight_kg marmite book 3 1.58 70.3 FALSE Life of Pi \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Creating and querying your own dataframe} Though used less, the \texttt{subset} command can be applied to extract information from selected columns only. For example: \begin{lstlisting} (small.survey2 = subset(small.survey, select = 3:4)) marmite book 1 TRUE Midnight in Chernobyl 2 FALSE To Kill a Mockingbird 3 FALSE Life of Pi \end{lstlisting} \pause And what if you wanted the height and favourite book of respondents who didn't like Marmite? \begin{lstlisting} subset(small.survey,marmite = = FALSE,select=c(1,4)) height_m book 2 1.75 To Kill a Mockingbird 3 1.58 Life of Pi \end{lstlisting} \end{frame} \subsubsection{The movies dataframe}\label{sec:moviesdataframe} \begin{frame}[fragile] \frametitle{The movies dataframe} The Movies dataset is available from within the \texttt{ncldata} package. \pause \vsp Recall that you can access this in the following way: \begin{lstlisting} install.packages("ncldata", repos="http://R-Forge.R-Project.org", type="source") library(ncldata) data(movies) \end{lstlisting} \pause \vspace{-0.5cm} Or: \begin{lstlisting} movies = read.table("https://www.mas.ncl.ac.uk/~nlf8/basicR/movies.txt", header=TRUE) \end{lstlisting} \vspace{-0.5cm} We will now explore the movies dataset in RStudio. \end{frame} \begin{frame} \frametitle{Example (One for later...)} \begin{tabular}{c l} \multirow{3}{*}{\includegraphics[width=1cm]{ironman.jpg}} &Load the \texttt{ncldata} package and read in the \texttt{movies} dataframe. \\ & Then answer the following questions. \end{tabular} \begin{itemize} \item [(a)] How many people submitted votes for the films ``\textit{10 Things I Hate About You}", ``\textit{Midnight in the Garden of Good and Evil}" and ``\textit{Zoolander}"? \item [(b)] How many films are classified as ``Action"? How many are classified as ``Romance"? How many ``Short" films are there in the database? \item [(c)] How many films were made between the Years 1970 and 2000 (inclusive)? \item [(d)] How many movies are longer than two hours in duration? \item [(e)] How many movies made since the year 2005 (not including 2005) are longer than three hours in duration? \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Example (One for later...)} \begin{tabular}{c l} \multirow{3}{*}{\includegraphics[width=1cm]{ironman.jpg}} & Load and explore the \texttt{music} dataset from the \texttt{ncldata} \\ &package by typing the following lines of code:\footnote{This will only work if you have correctly installed the \texttt{ncldata} package and read the package into your working RStudio space using the \texttt{library} command. Or, read in directly from my webpage, replacing \textcolor{blue}{\texttt{movies.txt}} with \textcolor{blue}{\texttt{music.txt}}.} \\ \end{tabular} \begin{lstlisting} data(music) #%* \textcolor{ForestGreen}{What's this dataset all about? Let's find out}*) ?music \end{lstlisting} \end{frame} \begin{frame} \frametitle{Example (One for later...)} \begin{itemize} \item [(a)] How many tracks do we have data for, and on how many variables? \item [(b)] How many tracks do we have from \textit{The Beatles} in this dataframe? What about \textit{Lady GaGa}? \item [(c)] How many tracks have a tempo of at least 130 beats per minute? \item [(d)] Of the tracks with a tempo of at least 130 beats per minute, how many have been released since the year 2000 (including the year 2000)? \item [(e)] What is the most popular key for tracks in this database? \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} R has very powerful built-in functions for summarising datasets both \textcolor{blue}{\bf graphically} and \textcolor{blue}{\bf numerically}. \vsp Most of these are fairly intuitive, so in this section we provide only a \textcolor{blue}{\bf brief review}. \vsp For example, we might want to explore the distribution of budgets in the \texttt{movies} dataframe. \vsp As we know, there are some missing values in the dataset, but most of the built-in functions for data analysis in R have an argument we can invoke to ignore these. \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{lstlisting} #%* \textcolor{ForestGreen}{To find the mean of the movie budgets but ignoring the missing values}*) mean(movies$Budget, na.rm = TRUE) [1] 41245023 #%* \textcolor{ForestGreen}{And the standard deviation}*) sd(movies$Budget, na.rm = TRUE) [1] 68612118 #%* \textcolor{ForestGreen}{How about a range of numerical summaries with a single function?}*) summary(movies$Budget, na.rm = TRUE) Min. 1st Qu. Median Mean 3rd Qu. Max. NA 0 9e+6 2.5e+7 4.125e+7 5e+7 1.9e+9 998 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} From the listing above we can make a few observations: \begin{itemize} \item The mean budget is larger than the median: $\sim$US\$ 41M, compared to US\$ 25M. This indicates a \textcolor{blue}{\bf positively skewed} dataset, as the mean is being `pulled up' by some very large budgets (of around USD\$1.9Bn?! Which movie is this?! Must be amazing...) \pause \item The standard deviation seems extremely large; again, this could be due to a small number of extremely large budgets... \pause \item ... or perhaps some unrealistically \textit{small} budgets (of US\$ 0?!) \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{lstlisting} #%* \textcolor{ForestGreen}{Let's find out the position of the most expensive movie}*) which(movies$Budget = = max(movies$Budget, na.rm=TRUE)) [1] 1211 #%* \textcolor{ForestGreen}{What's the \href{https://www.youtube.com/watch?v=PoPv-Ss6YbI}{\beamergotobutton{name}} of this movie?}*) movies$Title[1211] [1] "Enthiran" #%* \textcolor{ForestGreen}{Never heard of it!}*) #%* \textcolor{ForestGreen}{Which movies were made for free?!}*) movies$Title[which(movies$Budget==0)] [1] "Homefront" "Just Friends" "Knock Knock" [4] "Now You See Me" "Safe House" "Snitch" \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} What about some \textcolor{blue}{\bf graphical summaries}? \vsp By far the most commonly-used graphics are \textcolor{blue}{\bf histograms} and \textcolor{blue}{\bf boxplots}. \vsp Of course, there are so many more graphical options available in R; you should consult one of the books/texts mentioned at the start for a more comprehensive coverage. \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{lstlisting} #%* \textcolor{ForestGreen}{Partitions the plotting space: 1 row, 3 columns}*) par(mfrow = c(1,3)) #%* \textcolor{ForestGreen}{Now we'll try three different methods for producing histograms}*) #%* \textcolor{ForestGreen}{1. Sturges' rule; the default, tries to find an optimal number of bars}*) hist(movies$Budget, col="wheat",breaks="sturges", xlab="Budget (USD)", main="Sturges' rule") #%* \textcolor{ForestGreen}{2. Freedman-Diaconis's rule; tries to optimise the class width}*) hist(movies$Budget, col="gold",breaks="fd", xlab="Budget (USD)", main="FD rule") #%* \textcolor{ForestGreen}{3. Scott's rule; also tries to optimise the class width}*) hist(movies$Budget ,breaks="scott", xlab="Budget (USD)", main="Scott's rule") \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{center} %\begin{figure}[!p] %\centering \includegraphics[width=10.5cm, height=6cm]{../graphics/hist_2905.pdf} %\includegraphics[]{../graphics/boxplot_2905.pdf} %\caption{Histograms for the movie budgets, constructed using three different methods (top row); boxplots (bottom row), with outliers removed (right).} %\label{fig:hist} %\end{figure} \end{center} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{lstlisting} par(mfrow=c(1,2)) #%* \textcolor{ForestGreen}{Let's try a boxplot now}*) boxplot(movies$Budget, col="gold", main="Boxplot of movie budgets", horizontal=TRUE, xlab="Budget (USD)") #%* \textcolor{ForestGreen}{Outliers are spoiling things a bit. Let's remove them}*) boxplot(movies$Budget, col="gold", main="Boxplot of movie budgets", horizontal=TRUE, xlab="Budget (USD)", outline=FALSE) \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{center} %\begin{figure}[!p] %\centering %\includegraphics[]{../graphics/hist_2905.pdf} \includegraphics[width=10.5cm, height=6cm]{../graphics/boxplot_2905.pdf} %\caption{Histograms for the movie budgets, constructed using three different methods (top row); boxplots (bottom row), with outliers removed (right).} %\label{fig:hist} %\end{figure} \end{center} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} \begin{itemize} \item The histograms reveal just how positively skewed the movie budget data are \pause \item They also show how care should be taken with the method chosen to produce the histogram -- of course, this is a subjective choice, and some experimentation might be necessary \pause \item The boxplots also indicate just how skewed our dataset is -- to get a handle on the distribution of most of of our data, it might be preferable to produce the boxplot with the outliers removed \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Graphical and numerical summaries} %\begin{table}[!b] %\centering \begin{tabular}{l l l} \hline Command &Comment &Example \\ \hline \texttt{table} &Contingency table &\texttt{table(x)} \\ \texttt{barplot} &Generate a bar chart &\texttt{barplot(table(x))} \\ \texttt{hist} &Histogram &\texttt{hist(x)} \\ \texttt{plot} &Scatterplot &\texttt{plot(x, y)} \\ \texttt{points} &Add points to a plot &\texttt{points(x, y)} \\ \texttt{lines} &Add lines to a plot &\texttt{lines(x, y)} \\ \texttt{boxplot} &Box and whisker plot &\texttt{boxplot} \\ \hline \end{tabular} %\caption{Summary of R commands for commonly-used tables/graphs/graphical features} %\label{table:summaries} %\end{table} \end{frame} \begin{frame}[fragile] \frametitle{Example (One for later...)} \begin{tabular}{c l} \multirow{2}{*}{\includegraphics[width=1cm]{ironman.jpg} } & Load the \texttt{music} dataset into RStudio, just as you \\ &did earlier.\\ \end{tabular} \begin{itemize} \item [(a)] Find the mean and standard deviation for the duration of tracks in this dataset. \item [(b)] Produce histograms for track duration and danceability; put these histograms on the same panel of plots, use different colours for the histograms (you can choose whichever colours you like!) and make sure you provide suitable titles and labels on the axes. \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Example (One for later...)} \begin{itemize} \item [(c)] You are asked to investigate whether or not the duration of tracks has changed over time. \begin{itemize} \item [(i)] Use the \texttt{plot} command in R to produce a scatterplot of the year each track was released against track duration. \item [(ii)] Try out the following code, and comment on what you see in the plots in part (c): \begin{lstlisting} boxplot(music$duration~music$decade) \end{lstlisting} \vspace{-0.7cm} \end{itemize} \item [(d)] On average, how do energy scores for ``Hip Hop" tracks compare to ``Electronic" tracks? \end{itemize} \end{frame} \frame{ \frametitle{Functions, for loops and if statements} A very powerful aspect of R is that it is relatively easy to write your own functions. \vsp Functions take \textcolor{blue}{\bf inputs} and return a \textcolor{blue}{\bf single value}. \vsp Let's look at some simple functions. } \begin{frame}[fragile] \frametitle{Functions} \begin{lstlisting} Fun1 = function(x){ return(x*x) } \end{lstlisting} \pause The key elements here are: \begin{itemize} \item The word \textcolor{blue}{\bf function} \item The brackets () which enclose the \textcolor{blue}{\bf argument list} \item A sequence of statements inside braces $\{\}$ \item A \textcolor{blue}{\bf return} statement \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Functions} We \textcolor{blue}{\bf call} \texttt{Fun1} in the following way: \begin{lstlisting} Fun1(5) [1] 25 y = Fun1(10) y [1] 100 z = c(1, 2, 3, 4) Fun1(z) [1] 1 4 9 16 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Functions} Other variations include: \begin{lstlisting} #%* \textcolor{ForestGreen}{Default argument}*) Fun2 = function(x=1){ return(x*x) } Fun2() [1] 1 Fun2(4) [1] 16 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Functions} \begin{lstlisting} #%* \textcolor{ForestGreen}{Two arguments}*) Fun3 = function(x, y){ return(x*y) } Fun3(3,4) [1] 12 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Functions} When we call a function, R first looks for \textcolor{blue}{\bf local} variables and then \textcolor{blue}{\bf global} variables. \vsp For example, \texttt{Fun4} uses a global variable: \begin{lstlisting} blob = 5 Fun4 = function(){ return(blob) } Fun4() [1] 5 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Functions} However, \texttt{Fun5} uses a local variable: \begin{lstlisting} Fun5 = function(){ blob = 6 return(blob) } Fun5() [1] 6 blob [1] 5 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#4}} In the code below, what is the value of y? \begin{lstlisting} fun1 = function(x=3){ return(2*x) } y = fun1(2) \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 3 \item [2.] 4 \item [3.] 6 \item [4.] 10 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#4}} In the code below, what is the value of y? \begin{lstlisting} fun1 = function(x=3){ return(2*x) } y = fun1(2) \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 3 \item [2.] \textcolor{red}{\bf 4} \item [3.] 6 \item [4.] 10 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#5}} In the code below, what is the value of y? \begin{lstlisting} fun1 = function(x=3){ x = 2 return(2*x) } y = fun1() \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 3 \item [2.] 4 \item [3.] 6 \item [4.] 10 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#5}} In the code below, what is the value of y? \begin{lstlisting} fun1 = function(x=3){ x = 2 return(2*x) } y = fun1() \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 3 \item [2.] \textcolor{red}{\bf 4} \item [3.] 6 \item [4.] 10 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#6}} In the code below, what is the value of y? \begin{lstlisting} x = 1 fun1 = function(x){ x = 2 return(x) } y = fun1(3) \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 1 \item [2.] 2 \item [3.] 3 \item [4.] 4 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#6}} In the code below, what is the value of y? \begin{lstlisting} x = 1 fun1 = function(x){ x = 2 return(x) } y = fun1(3) \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 1 \item [2.] \textcolor{red}{\bf 2} \item [3.] 3 \item [4.] 4 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#7}} In the code below, what is the value of y? \begin{lstlisting} x = 1 fun1 = function(x = 3){ x = 2 return(x) } y = fun1() \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 1 \item [2.] 2 \item [3.] 3 \item [4.] 4 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#7}} In the code below, what is the value of y? \begin{lstlisting} x = 1 fun1 = function(x = 3){ x = 2 return(x) } y = fun1() \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} \begin{itemize} \item [1.] 1 \item [2.] \textcolor{red}{\bf 2} \item [3.] 3 \item [4.] 4 \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#8}} In the code below, what is the value of y? \begin{lstlisting} x = 1 fun1 = function(x = 3){ return(2*x) } x = fun1(x) y = fun1(x) \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} {\huge \textcolor{white}{\bf 4}} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#8}} In the code below, what is the value of y? \begin{lstlisting} x = 1 fun1 = function(x = 3){ return(2*x) } x = fun1(x) y = fun1(x) \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} {\huge \textcolor{red}{\bf 4}} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: If statements} Conditional statements are features of a programming language which perform different computations or actions depending on whether a condition evaluates \texttt{TRUE} or \texttt{FALSE}. \vsp The basic structure of an \textcolor{blue}{\bf if statement} is: \begin{lstlisting} if(expr) { ##do something } \end{lstlisting} where \texttt{expr} is evaulated to be either \texttt{TRUE} or \texttt{FALSE}. \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: If statements} For example: \begin{lstlisting} x = 5 y = 5 if(x < 5){ y = 0 } y [1] 5 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: If statements} Another example: \begin{lstlisting} x = 5 y = 5 if(x > 0){ y = 0 } y [1] 0 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: If else statements} We can link together a number of if statements: \pause \begin{lstlisting} x = 0 if(x > 0){ cat("x is greater than zero") } else if(x < 0){ cat("x is less than zero") } else{ cat("x must be zero!") } x must be zero! \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: If statements inside functions} We can also use if statements inside functions. For example: \begin{lstlisting} IsNegative = function(value){ I = FALSE if(value < 0){ I = TRUE } return(I) } IsNegative(1) [1] FALSE IsNegative(-5.6) [1] TRUE \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#9}} In the code below, what is the value of z? \begin{lstlisting} x = 5 y = "male" if(x > 6 & y == "Female"){ z = 0 } else if(x < 5 & y =="Male"){ z = 1 } else{ z = 2} \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} {\huge \textcolor{white}{\bf 2}} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#9}} In the code below, what is the value of z? \begin{lstlisting} x = 5 y = "male" if(x > 6 & y == "Female"){ z = 0 } else if(x < 5 & y =="Male"){ z = 1 } else{ z = 2} \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} {\huge \textcolor{red}{\bf 2}} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: For loops} \begin{itemize} \item At times we would like to perform some operation on a vector or a data frame \item Often R has built-in functions that will do this for you, e.g. \texttt{mean}, \texttt{sd}; other times we have to write our own functions \item For example, suppose we want to calculate $\sum_{i=1}^{10} i^2 \;$ \item In R we can use a \texttt{for} loop: \pause \end{itemize} \begin{lstlisting} x = 0 for(i in 1:10){ x = x + i^2 } \end{lstlisting} \pause \begin{lstlisting} x [1] 385 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: For loops} Example: Find \[ \sum_{j=-5}^{-1} e^{j}/j^{2}. \] \pause \begin{lstlisting} total = 0 for(j in -5:-1){ total = total + exp(j)/j^2 } \end{lstlisting} \pause \begin{lstlisting} total [1] 0.4086594 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Conditionals: For loops} Example: Find \[ \sum e^{k}/k^{2}, \quad k=3, 6, 9, \ldots, 21. \] \pause \begin{lstlisting} total = 0 for(i in 1:7){ k = i*3 total = total + exp(k)/k^2 } total [1] 3208939 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#10}} What is the value of x? \begin{lstlisting} total = 0 for(blob in 1:10){ if(blob>9){ total = total + 1 } } x = total \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} {\huge \textcolor{white}{\bf 1}} \end{frame} \begin{frame}[fragile] \frametitle{\textcolor{purple}{Quiz \#10}} What is the value of x? \begin{lstlisting} total = 0 for(blob in 1:10){ if(blob>9){ total = total + 1 } } x = total \end{lstlisting} \includegraphics[width=1cm]{ironman.jpg} {\huge \textcolor{red}{\bf 1}} \end{frame} \end{document} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% OLD STUFF \begin{frame}\frametitle{MAS2602 Arrangements} \begin{itemize} \item Statistics classes run in teaching weeks 5--9 \item Lecturer is Dr. Lee Fawcett (\textcolor{blue}{\texttt{lee.fawcett@ncl.ac.uk}}) \item Schedule: 4 lectures, 4 practicals, 1 revision class, 1 class test \item Office hours: Tuesdays 3-4; room 2.07 Herschel \\ (also MAS3902: Tuesdays 2-3 and Thursdays 3-4) \end{itemize} \end{frame} \begin{frame} \frametitle{Schedule} %\vspace{-1.1cm} {\small \begin{center} \begin{tabular}{llll} Week 5 &Mon 28 Oct &11--12 &Lecture (Herschel LT1) \\ &Fri 1 Nov &10--12 &Practical (Herschel PC) \\ & & & \\ Week 6 &Mon 4 Nov &11--12 &Lecture (Herschel LT1) \\ &Wed 6 Nov &11--1 &Practical (Herschel PC) \\ & & & \\ Week 7 &Mon 11 Nov &11--12 &Lecture (Herschel LT1) \\ &Thu 14 Nov &11--1 &Practical (Herschel PC) \\ & & & \\ Week 8 &Mon 18 Nov &11--12 &Lecture (Herschel LT1) \\ &Thu 21 Nov &11--1 &Practical (Herschel PC) \\ &\textcolor{red}{\it Thu 21 Nov} &\textcolor{red}{\it 3pm} &\textcolor{red}{\it Assignment due} \\ &Fri 22 Nov &1--2 &Revision (Herschel LT2) \\ & & & \\ Week 9 &\textcolor{red}{\it Tue 26 Nov} &\textcolor{red}{\it 9--10} &\textcolor{red}{\it Test (Herschel PC)} \\ \end{tabular} \end{center} } \end{frame} \begin{frame} \frametitle{Assessment} \textbf{Assignment (a.k.a.\ ``mini-project'')} \begin{itemize} \item Due \textcolor{red}{\bf Thursday 21 November, 3pm} \item Worth $10\%$ of the module marks \end{itemize} \textbf{Class test} \begin{itemize} \item \textcolor{red}{\bf Tuesday 26 November, 9am} -- one hour long \item Worth $30\%$ of module marks \item Open-book test \item You will write one or two short R programs during the test \end{itemize} \end{frame} \begin{frame} \frametitle{Help and Support} Help available: \begin{itemize} \item Office hours \item Demonstrators in practical sessions \item Books -- see recommendations in booklet \item Email Lee, or just pop in! \item Blackboard and dedicated webpage \end{itemize} \end{frame} \begin{frame} \frametitle{Late work policy} \begin{center} In this module, deadline extensions can be requested for the final project (by means of submitting a PEC form), and work submitted within 7 days of the deadline without good reason will be marked for reduced credit. This module also contains tests worth more than $10\%$ for which rescheduling can be requested (by means of submitting a PEC form). There are mini-projects (worth $10\%$ each) for which it is not possible to extend deadlines and for which no late work can be accepted. For details of the policy (including procedures in the event of illness etc.) please look at the School web site: http://www.ncl.ac.uk/maths/students/teaching/homework/ For problems with deadlines, speak to your personal tutor and prepare a PEC form \end{center} \end{frame} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{frame} \begin{Large} \begin{center} {\bf Lecture 1: Introduction and Simulation of Random Variables} \end{center} \end{Large} \end{frame} \begin{frame}\frametitle{Introduction} In this part of the module we will do statistics with R: \begin{itemize} \item R is the foremost tool in modern computational statistics \item Using R teaches general concepts in programming \item It can be used to illustrate mathematical ideas in probability and statistics \end{itemize} \pause \textbf{Today's lecture:} simulating random variables \begin{enumerate} \item Simulating random variables seen in MAS1604 \item Using this to simulate more complicated probability models \end{enumerate} \pause \textbf{Friday's practical:} \begin{enumerate} \item Revision of R from MAS1608 \item Putting today's material into practice \end{enumerate} \end{frame} \begin{frame}\frametitle{The binomial distribution} If $X\sim Bin(n,p)$ then $X$ has PMF\ (probability mass function) given by \textcolor{BrickRed}{ \begin{equation*} p_X(k) = \binom{n}{k}p^k(1-p)^{n-k},\quad\text{for\ }k=0,1,\ldots,n. \end{equation*}} There is no closed formula for the CDF (cumulative distribution function). \pause R commands: \begin{itemize} \item \texttt{dbinom} -- calculate PMF \item \texttt{pbinom} -- calculate CDF \item \texttt{rbinom} -- generate random sample \end{itemize} \end{frame} \begin{frame}[fragile]\frametitle{R commands for the binomial distribution} \begin{lstlisting} dbinom(5, 10, 0.7) #%* \textcolor{ForestGreen}{$\pr{X=5}, X\sim Bin(10,0.7)$} *) pbinom(4, 10, 0.7) #%* \textcolor{ForestGreen}{$\pr{X\leq 4}, X\sim Bin(10,0.7)$} *) rbinom(50, 10, 0.7) #%* \textcolor{ForestGreen}{Sample a $Bin(10,0.7)$ distribution 50 times} *) \end{lstlisting} \end{frame} \begin{frame}[fragile]\frametitle{Creating a bar plot} \begin{lstlisting} x = rbinom(50, 10, 0.7) plot(table(x), xlim = c(0,10), xaxt = 'n', xlab = 'x', ylab = 'frequency') axis(1, at = seq(0, 10)) \end{lstlisting} \vspace{-1.5cm} \begin{center} \includegraphics[width=8cm]{../notes/fig1-1} \end{center} \end{frame} \begin{frame}\frametitle{The geometric distribution} If $Y\sim Geom(p)$ then $Y$ has PMF and CDF given by \textcolor{BrickRed}{ \begin{align*} p_Y(k) &= (1-p)^{k-1}p,\\ F_Y(k) &= 1-(1-p)^k,\quad\text{for\ }k=1,2,\ldots. \end{align*}} Note that $Y$ takes values in $1,2,3,\ldots$: \begin{itemize} \item R uses a slightly different definition \item We use the definition that $Y$ is the number of Bernoulli trials with up to and including first success \item R counts number of trials up to, but \hi{not} including the first success, so in R geometric random variables take values $0,1,2,\ldots$. \end{itemize} \end{frame} \begin{frame}[fragile]\frametitle{R commands for the geometric distribution} Adjust the arguments to account for different definition: \begin{lstlisting} dgeom(4, 0.2) #%* \textcolor{ForestGreen}{$\pr{Y=5}, Y\sim Geom(0.2)$} *) pgeom(2, 0.2) #%* \textcolor{ForestGreen}{$\pr{Y\leq 3}, Y\sim Geom(0.2)$} *) 1 + rgeom(100, 0.2) #%* \textcolor{ForestGreen}{Sample a $Geom(0.2)$ distribution 100 times} *) \end{lstlisting} Here's a function to replace \texttt{dgeom} with our definition of the geometric distribution: \begin{lstlisting} mydgeom = function(x, p) { dgeom(x-1, p)} \end{lstlisting} \end{frame} \begin{frame}\frametitle{The Poisson distribution} If $Z\sim Po(\lambda)$ then $Z$ has PMF given by \textcolor{BrickRed}{ \begin{equation*} p_Z(k) = \frac{\lambda^k}{k!}e^{-\lambda},\quad\text{for\ }k=0,1,2,\ldots. \end{equation*}} There is no closed formula for the CDF (cumulative distribution function). R commands: \begin{itemize} \item \texttt{dpois} -- calculate PMF \item \texttt{ppois} -- calculate CDF \item \texttt{rpois} -- generate random sample \end{itemize} \end{frame} \begin{frame}[fragile]\frametitle{R commands for the Poisson distribution} \begin{lstlisting} dpois(5, 3.5) #%* \textcolor{ForestGreen}{$\pr{Z=5}, Z\sim Po(3.5)$} *) ppois(2, 3.5) #%* \textcolor{ForestGreen}{$\pr{Z\leq 2}, Z\sim Po(3.5)$} *) rpois(100, 3.5) #%* \textcolor{ForestGreen}{Sample a $Po(3.5)$ distribution $100$ times} *) \end{lstlisting} \vspace{-0.5cm} \end{frame} \begin{frame}\frametitle{Summary} \begin{center} \begin{tabular}{c|ccc} Distribution & Binomial & Poisson & Geometric\ $\skull$\\ \hline PMF & \texttt{dbinom(...)} & \texttt{dpois(...)} & \texttt{dgeom(...)} \\ CDF & \texttt{pbinom(...)} & \texttt{ppois(...)} & \texttt{pgeom(...)} \\ sample & \texttt{rbinom(...)} & \texttt{rpois(...)} & \texttt{rgeom(...)} \\ \end{tabular} \end{center} \end{frame} \begin{frame}\frametitle{Continuous random variables} R has functions for the uniform, exponential and normal distributions: \begin{center} \begin{tabular}{c|ccc} Distribution & Uniform & Exponential & Normal\\ \hline PDF & \texttt{dunif(...)} & \texttt{dexp(...)} & \texttt{dnorm(...)} \\ CDF & \texttt{punif(...)} & \texttt{pexp(...)} & \texttt{pnorm(...)} \\ quantile & \texttt{qunif(...)} & \texttt{qexp(...)} & \texttt{qnorm(...)} \\ sample & \texttt{runif(...)} & \texttt{rexp(...)} & \texttt{rnorm(...)} \\ \end{tabular} \end{center} \end{frame} \begin{frame}[fragile]\frametitle{Continuous random variables -- R examples} \begin{lstlisting} dunif(3, 2, 5) #%* \textcolor{ForestGreen}{$f_{X}(3), X\sim U(2, 5)$} *) pexp(1.5, 5) #%* \textcolor{ForestGreen}{$F_{Y}(1.5), Y\sim Exp(5)$} *) rnorm(10, 0, 2) #%* \textcolor{ForestGreen}{Sample a $N(0, 2^2)$ distribution 10 times} *) \end{lstlisting} For the standard uniform ($U(0,1)$) and standard normal ($N(0,1)$) distributions you don't need to provide the parameters $a=0$, $b=1$ and $\mu=0$, $\sigma=1$ respectively. For example: \begin{lstlisting} runif(20) #%* \textcolor{ForestGreen}{Samples a $U(0,1)$ distribution 20 times} *) pnorm(1.96) #%* \textcolor{ForestGreen}{$\pr{Z<1.96}, Z\sim N(0,1)$} *) [1] 0.9750021 \end{lstlisting} \end{frame} \begin{frame}[fragile]\frametitle{Quantiles} The quantile functions \texttt{qunif, qexp, qnorm} solve equations like \begin{equation*} F_X(\alpha) = p \end{equation*} for $\alpha$ given a probability $p$. For example: \begin{lstlisting} qnorm(0.9750021) [1] 1.96 \end{lstlisting} \end{frame} \begin{frame}\frametitle{Quantiles} \begin{example} Suppose annual maximum wave heights observed off the coast at a flood-prone town are assumed Normally distributed, with mean 2 metres and standard deviation 0.5 metres. \begin{itemize} \item [(a)] Write down the R command to find the probability that the largest wave height next year will exceed 3.25 metres. \item [(b)] Write down the R command to estimate the height of a new sea wall such that we might expect the town to be flooded, on average, once per century. Why might our modelling assumption be invalid? \end{itemize} \end{example} \end{frame} \frame{ \frametitle{Example 1.1: Solution (a)} \includegraphics[width=\textwidth, height=\textheight]{seasurge} } \begin{frame}[fragile] \frametitle{Example 1.1: Solution (a)} Could work out the old-fashioned way: \pause \begin{eqnarray*} \pr{X>3.25} &=&\pr{Z>\frac{3.25-2}{0.5}} \\ & & \\ &=& \pr{Z>2.5} \\ & & \\ &=& 1-\pr{Z \leq 2.5} = 1-0.994 = 0.006. \end{eqnarray*} \pause Or could just use R: \begin{lstlisting} 1-pnorm(3.25, 2, 0.5) [1] 0.006209665 \end{lstlisting} \end{frame} \begin{frame}[fragile] \frametitle{Example 1.1: Solution (b)} \includegraphics[width=\textwidth, height=\textheight]{seasurge2} \end{frame} \begin{frame}[fragile] \frametitle{Example 1.1: Solution (b)} Similarly: \begin{lstlisting} qnorm(0.99, 2, 0.5) [1] 3.163174 \end{lstlisting} \end{frame} \begin{frame}\frametitle{A more advanced model} \begin{itemize} \item Number of arrivals per day at an IT help-desk is modelled using a Poisson distribution \item Mean of the Poisson distribution might vary from day-to-day \item Suppose the number of arrivals $X\sim Po(\Lambda)$ \item $\Lambda$ is itself a random variable, with $\Lambda\sim Exp(c)$ for a constant $c=0.05$ \end{itemize} What can we say about the distribution of $X$? \begin{enumerate} \item What are the expectation and variance of $X$? \item What is $\pr{X>30}$? \end{enumerate} \end{frame} \begin{frame}[fragile]\frametitle{Arrival model -- R code} \begin{lstlisting} arrivals = function(n, c = 0.05) { x = vector(mode = 'numeric',length = n) for (i in 1:n) { lambda = rexp(1, rate = c) x[i] = rpois(1, lambda) } return(x) } \end{lstlisting} \end{frame} \begin{frame}[fragile]\frametitle{Arrival model -- R code} Bar plot: \begin{lstlisting} x = arrivals(500) plot(table(x), xlim = c(0, 50), xaxt = 'n', xlab = 'x', ylab = 'frequency') axis(1, at = seq(0, 50, 5)) \end{lstlisting} \vspace{-1.5cm} \begin{center} \includegraphics[width=8cm]{../notes/fig1-2} \end{center} \end{frame} \begin{frame}\frametitle{Using simulated samples} \begin{itemize} \item The \hi{sample mean} is an approximation to the distribution mean \item The same applies for the sample variance \item To calculate $\pr{X> 30}$ we count the proportion of times this occurs in the sample \end{itemize} We expect the approximation to improve \hi{as we increase the sample size} \end{frame} \begin{frame}[fragile]\frametitle{Using simulated samples} \begin{lstlisting} mean(x) [1] 19.618 var(x) [1] 382.1764 sum(x>30)/500 [1] 0.202 \end{lstlisting} \end{frame} \begin{frame}\frametitle{Mixtures of normal distributions} Suppose \begin{itemize} \item $\mu_1,\mu_2,\sigma_1>0,\sigma_2>0$ are fixed constants, \item $w_1,w_2$ are positive constants with $w_1+w_2=1$. \end{itemize} Consider the following function: \begin{equation*} f(x) = w_1f_1(x)+w_2f_2(x) \end{equation*} where $f_1$ and $f_2$ are the density functions for $Z_1\sim N(\mu_1,\sigma_1^2)$ and $Z_2\sim N(\mu_2,\sigma_2^2)$ respectively. Check that $f(x)$ represents a valid probability density function. \end{frame} \begin{frame}\frametitle{Mixtures of normal distributions} \begin{answers} First note that $f(x) \geq 0$ everywhere. Also \begin{eqnarray*} \int_{-\infty}^{\infty} f(x) dx &=& w_{1}\int_{-\infty}^{\infty}f_{1}(x) dx + w_{2} \int_{-\infty}^{\infty}f_{2}(x)dx \pause \\ & & \\ &=& w_1 + w_2 = 1. \end{eqnarray*} \end{answers} \end{frame} \begin{frame}\frametitle{Mixtures of normal distributions} If a random variable $X$ has PDF corresponding to $f(x)$ we say it is a \hi{mixture} of normal distributions $N(\mu_1,\sigma_1^2)$ and $N(\mu_2,\sigma_2^2)$ with weights $w_1,w_2$. Note that $X$ is \hi{not} the sum of two normal random variables i.e. \begin{equation*} X \neq w_1Z_1+w_2Z_2\quad\text{where\ }Z_i\sim N(\mu_i,\sigma_i^2),\ i=1,2. \end{equation*} \end{frame} \begin{frame}\frametitle{Example} For example: suppose $\mu_1=3, \sigma_1=1$ and $\mu_2=6, \sigma_2=2$ with $w_1=w_2=1/2$. \vspace{-1.7cm} \begin{center} \includegraphics[width=9.5cm]{../notes/fig1-3} \end{center} \end{frame} \begin{frame}\frametitle{Sampling from a mixture distribution} \begin{enumerate} \item Sample a random variable $J\in\{1,2\}$ such that $\pr{J=1}=w_1$ and $\pr{J=2}=w_2$. \item The random variable $J$ tells you which \hi{component} of the mixture to sample $X$ from. \item If $J=1$ then sample $X$ from $N(\mu_1,\sigma_1^2)$, but if $J=2$ then sample $X$ from $N(\mu_2,\sigma_2^2)$. \end{enumerate} \end{frame} \begin{frame}[fragile]\frametitle{Sampling from a mixture distribution} \begin{lstlisting} normal.mixture = function(n, mu1, sig1, w1, mu2, sig2, w2) { p = c(w1, w2) x = vector(mode = 'numeric', length = n) for (i in 1:n) { j = sample(c(1, 2), 1, prob = p) if (j==1) { x[i] = rnorm(1, mu1, sig1) } else { x[i] = rnorm(1, mu2, sig2) } } return(x) } \end{lstlisting} \end{frame} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{frame} \begin{Large} \begin{center} {\bf Lecture 2: Inverse CDF method} \end{center} \end{Large} \end{frame} \begin{frame}\frametitle{Today's lecture} In this lecture we look at the \textcolor{purple}{\bf inverse CDF method} for simulating from continuous random variables. \end{frame} \begin{frame}\frametitle{The inverse CDF method} This is a method for simulating \hi{univariate continuous} random variables \begin{itemize} \item Let $U\sim U(0,1)$ \item Suppose $F(x)$ is a well-defined CDF which is invertible \item Then the random variable $X=F^{-1}(U)$ has CDF $F(x)$ \end{itemize} \end{frame} \begin{frame} \frametitle{The inverse CDF method} \begin{answers} \[ X = F^{-1}(U). \pause \] So \begin{eqnarray*} \pr{X \leq x} &=& \pr{F^{-1}(U) \leq x} \\ \pause & & \\ &=&\pr{F(F^{-1}(U)) \leq F(x)} \\ \pause & & \\ & = &\pr{U \leq F(x)}\\ \pause & &\\ &=&F(x), \end{eqnarray*} since $0 \leq F(x) \leq 1$ (see sketch on board). \end{answers} \end{frame} \begin{frame}\frametitle{The inverse CDF method} \begin{example} Suppose we obtain two observations from a $U(0,1)$ distribution: $0.1$ and $0.85$. Use these values to obtain two observations from $X\sim Exp(2)$. \end{example} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} \noindent If $X \sim Exp(2)$ then $F_{X}(x) = 1-e^{-2x}$. \pause \vspace{0.3cm} If $y = 1-e^{-2x}$, then \begin{eqnarray*} e^{-2x} &=& 1-y \\ & & \\ -2x &=& \log (1-y) \\ & & \\ x &=& -\frac{1}{2} \log(1-y). \end{eqnarray*} \pause Plug in $y=0.1$ and $y=0.85$; gives realized values of $0.0527$ and $0.949$. \end{answers} \end{frame} \begin{frame}\frametitle{Example} \begin{itemize} \item Patients arriving at a doctor's surgery have appointment times $15$ minutes apart. \item The random variable $L$ measures how late a patient is for his/her appointment, with negative values denoting early arrivals. \end{itemize} The PDF of $L$ is \begin{equation*} f_L(\ell) = \begin{cases} 0.0025(\ell+20),\ & -20\leq \ell \leq 0,\\ 0.05e^{-\ell/10},& \ell>0. \end{cases} \end{equation*} \end{frame} \begin{frame}\frametitle{Questions} \begin{enumerate}[a.] \item Sketch $f_L(\ell)$. \item Show that the CDF of $L$ is \begin{equation*} F_L(\ell) = \begin{cases} 0,&\ell<-20,\\ 0.00125(\ell+20)^2, & -20\leq \ell \leq 0,\\ 1-\frac{e^{-\ell/10}}{2}, & \ell>0. \end{cases} \end{equation*} \item Suppose that $0.205$ and $0.713$ are two realized values of $U\sim U(0,1)$. Show that the corresponding realized values of $L$ are $-7.19$ and $5.55$ under the inverse CDF simulation method. \end{enumerate} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} \textbf{First, the PDF:} \vspace{-1.4cm} \begin{center} \includegraphics[width=8.25cm]{Ex24.pdf} \end{center} \end{answers} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} \noindent \textbf{Now to find the CDF:} \pause $L$ lies between $-20$ and $+\infty$, so $F_{L}(\ell) = 0$ when $\ell < -20$. \pause \vspace{0.3cm} \noindent For $-20 \leq \ell \leq 0$: \pause \begin{eqnarray*} F_{L}(\ell) &=& \int_{-20}^{\ell} f_{L}(s) ds \pause \\ & & \\ &=& \int_{-20}^{\ell} 0.0025(s+20) ds \\ \pause & & \\ &=&\left[\frac{0.0025}{2}(s+20)^{2}\right]^{\ell}_{-20} \pause\\ & & \\ &=& 0.00125(\ell+20)^{2} \quad \textsf{for $-20 \leq \ell \leq 0$}. \end{eqnarray*} Also, $F_{L}(0)=0.5$. \end{answers} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} Now \begin{eqnarray*} F_{L}(\ell) &=& \int_{-20}^{0}f_{L}(s) ds + \int_{0}^{\ell}f_{L}(s) ds \quad \textsf{for $\ell>0$} \pause \\ && \\ &=& 0.5 + \int_{0}^{\ell}0.05e^{-s/10} ds \pause \\ && \\ &=& 0.5 + \left[0.05 \times -10 e^{-s/10}\right]_{0}^{\ell} \pause \\ && \\ &=&0.5 +0.5 - 0.5 e^{-\ell/10} \pause \\ & & \\ &=& 1- \frac{e^{-\ell/10}}{2}, \quad \textsf{as required.} \end{eqnarray*} \end{answers} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} \noindent \textbf{Now for some realized values:} \pause \vspace{0.3cm} $u=0.205 <0.5$, so realized value $\ell$ lies between $-20$ and $0$. Thus \pause \[ F_{L}(\ell) = 0.00125(\ell+20)^2. \] \pause Find the inverse, and evaluate at $u=0.205$: \pause \begin{eqnarray*} y &=&0.00125(\ell+20)^2, \quad \textsf{so} \pause\\ \ell &=&\pm\sqrt{\frac{y}{0.00125}}-20, \end{eqnarray*}\pause giving \[ \ell = \sqrt{\frac{0.205}{0.00125}}-20 = -7.2. \] \end{answers} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} Similarly, $u=0.713 >0.5$ so realized value $>0$. Thus\pause \[ F_{L}(\ell) = 1- \frac{e^{-\ell/10}}{2}. \]\pause Find the inverse, and evaluate at $u=0.713$: \pause \begin{eqnarray*} y &=& 1-\frac{e^{-\ell/10}}{2} \quad \textsf{so} \pause \\ \ell &=& -10\log(2-2y), \end{eqnarray*} giving \[ \ell = -10\log(2-2\times 0.713) = 5.55. \] \end{answers} \end{frame} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{frame} \begin{Large} \begin{center} {\bf Lecture 3: Transformation of \\ Random Variables} \end{center} \end{Large} \end{frame} \begin{frame} \frametitle{Recap Quiz} \begin{itemize} \item [1.] Write down R code that simulates 100 observations from the $U(0,1)$ distribution and stores the output in the vector \texttt{U}. \item [2.] Suppose $X$ has CDF \[ F_{X}(x) = 1- \frac{1}{x}, \quad x\geq 1. \] Find $F^{-1}_{X}(x)$. \item [3. ] Write an R function called \texttt{inv.cdf} to return the inverse of $F_{X}(x)$. \item [4. ] Write down a single line of R code to generate 100 realizations of the random variable $X$. \item [5.] Write down a single line of R code to produce a histogram of your realizations of $X$, coloured yellow and with a suitable title. \end{itemize} \end{frame} \begin{frame}\frametitle{Transforming a normal random variable} \begin{itemize} \item Consider how we might model a \hi{continuous positive} random quantity e.g.\ blood pressure or height \item The Normal distribution is not suitable! \end{itemize} One way to obtain a positive random variable is to define \begin{equation*} Y = \exp(X),\quad\text{where\ }X\sim N(\mu,\sigma^2). \end{equation*} Then $Y$ is a random variable which is a \hi{transformation} of $X$. \end{frame} \begin{frame}[fragile]\frametitle{Simulation to explore transformation} \begin{lstlisting} y = exp(rnorm(250, mean = 1.5, sd = 0.3)) hist(y, xlab = 'Y', ylab = 'frequency') \end{lstlisting} \vspace{-1cm} \begin{center} \includegraphics[width=8cm]{../notes/fig3-1} \end{center} \end{frame} \begin{frame}\frametitle{Calculating PDFs} \begin{itemize} \item This looks like a useful distribution -- what is its PDF? \item $Y=\exp(X)$ so $X = \log(Y)$ \item So $f_Y(y)=f_X(\log(y))$ is a good guess, but is \hi{NOT} the right answer! \end{itemize} \vspace{1cm} \begin{example} If $X\sim U(0,1)$, what is the PDF of $Y=X^3$? \end{example} \end{frame} \begin{frame} \frametitle{Solution} \begin{answers} \begin{eqnarray*} F_{Y}(y) &=&\pr{Y \leq y} \pause\\ &=&\pr{X^{3} \leq y} \pause\\ &=& \pr{X \leq y^{1/3}} \\ &=& F_{X}(y^{1/3}). \pause \end{eqnarray*} Therefore \begin{eqnarray*} f_{Y}(y) &=& \frac{d}{dy} F_{X}(y^{1/3}) \pause \\ &=&f_{X}(y^{1/3}) \times \frac{1}{3}y^{-2/3} \quad \textsf{by chain rule} \pause \\ &=& \frac{1}{3}y^{-2/3}, \quad 0\leq y \leq 1 \textsf{ (and zero otherwise)}. \end{eqnarray*} \end{answers} \end{frame} \begin{frame}[fragile]\frametitle{Checking via simulation} To confirm this result, we can simulate $Y$: \begin{lstlisting} y = (runif(500))^3 hist(y, xlab = 'Y', ylab = 'frequency') \end{lstlisting} \vspace{-1cm} \begin{center} \includegraphics[width=8cm]{../notes/fig3-2} \end{center} \end{frame} \begin{frame}\frametitle{Further examples} \begin{example} If $X\sim N(0,1)$, what is the PDF of $Y=X^2$? \end{example} \pause \begin{answers} \[ F_{Y}(y) =\pr{Y \leq y} =\pr{X \leq y^{1/2}} = F_{X}(y^{1/2}). \] \pause Then \begin{eqnarray*} f_{Y}(y) &=& \frac{d}{dy} F_{X}(y^{1/2}) \pause \\ &=&f_{X}(y^{1/2}) \times \frac{1}{2}y^{-1/2} \quad \textsf{by chain rule} \pause \\ &=&\frac{1}{\sqrt{2\pi}}e^{-y/2} \times \frac{1}{2}y^{-1/2} \pause \\ &=&\frac{1}{2 \sqrt{2 \pi y}}e^{-y/2} \quad \textsf{for $y \geq 0$, zero otherwise.} \end{eqnarray*} \end{answers} \end{frame} \begin{frame}\frametitle{Further examples} \begin{example} If $X\sim N(\mu,\sigma^2)$, what is the PDF of $Y=\exp(X)$? \\ This is the example used at the start of the lecture. We say $Y$ has the \hi{log-normal} distribution and write $Y\sim LN(\mu,\sigma^2)$. \end{example} \pause \begin{answers} \begin{eqnarray*} F_{Y}(y) &=& F_{X}(\log y) \quad \textsf{by usual argument with CDF} \pause \\ f_{Y}(y) &=&\frac{d}{dy}F_{X}(\log y) \pause \\ &=&f_{X}(\log y) \times \frac{1}{y} \quad \textsf{for $y \geq 0$} \pause \\ &=&\frac{1}{\sigma y \sqrt{2 \pi}}\exp\left\{-\frac{1}{2}\left(\frac{(\log y - \mu)^{2}}{\sigma^{2}}\right)\right\}. \end{eqnarray*} \end{answers} \end{frame} \begin{frame}\frametitle{The univariate transformation method} If $X$ has pdf $f_X(x)$ and $h(x)$ is a monotonically increasing or monotonically decreasing function then \begin{itemize} \item $h$ is invertible and \item the PDF of $Y=h(X)$ is \begin{equation*} f_Y(y) = f_X(h^{-1}(y))\times\left| \frac{d\left( h^{-1}(y) \right)}{dy} \right| \end{equation*} when $y$ is in the range of $h(x)$ and zero otherwise. \end{itemize} This result allows direct calculation of $f_Y(y)$ and is called the \hi{univariate transformation method}. \begin{center} \textcolor{BrickRed}{\bf See whiteboard for Proof.} \end{center} \end{frame} \begin{frame}\frametitle{Some Examples} \begin{example} When $X\sim U(0,1)$, use the transformation method to find the PDF of $Y=X^3$. \end{example} \pause \begin{answers} \[ h(x) =x^{3}, \quad h^{-1}(y) = y^{1/3}. \] \pause \begin{eqnarray*} f_{Y}(y) &=& 1 \times \frac{d}{dy}y^{1/3} \pause \\ &=& \frac{1}{3}y^{-2/3}, \quad 0\leq y \leq 1. \end{eqnarray*} \end{answers} \end{frame} \begin{frame}\frametitle{Some examples} \begin{example} When $X\sim U(-1,1)$ and $Y=\exp(X)$, use the transformation method to find the PDF of $Y$. \end{example} \pause \begin{answers} \[ h(x) = e^{x}, \quad h^{-1}(y) = \log (y). \]\pause \begin{eqnarray*} f_{Y}(y) &=& \frac{1}{2} \times \abs[\Big]{\frac{d}{dy} \log (y)} \pause \\ &=& \frac{1}{2|y|} = \frac{1}{2y}, \quad e^{-1} \leq y \leq e^{1}, \textsf{ zero otherwise.} \end{eqnarray*} \end{answers} \end{frame} \begin{frame}\frametitle{Some examples} \begin{example} Use the transformation method to find the PDF of the log-normal distribution. \end{example} \pause \begin{answers} \[ h(x) = e^{x}, \quad h^{-1}(y) = \log (y). \] \pause \begin{eqnarray*} f_{Y}(y) &=& f_{X}(\log (y)) \times \frac{1}{y} \qquad \textsf{for $y>0$} \pause \\ &=& \frac{1}{y} \times \frac{1}{\sigma \sqrt{2 \pi}} \exp\left[-\frac{1}{2} \left(\frac{\log (y) - \mu}{\sigma}\right)^{2}\right]. \end{eqnarray*} \end{answers} \end{frame} \begin{frame}\frametitle{Cauchy example} \begin{example} Suppose $X\sim U(-\frac{\pi}{2}, \frac{\pi}{2})$ and $Y=\tan(X)$. \\ Use simulation to study the distribution of $Y$ and explicitly compute the PDF. \end{example} \end{frame} \begin{frame}[fragile]\frametitle{Solution} \begin{lstlisting} y = tan(runif(100, -pi/2, pi/2)) hist(y, xlab = 'y', ylab = 'frequency', main = '') \end{lstlisting} \vspace{-1cm} \begin{center} \includegraphics[width=8cm]{../notes/fig3-3} \end{center} \end{frame} \begin{frame}\frametitle{Solution} \begin{answers} This distribution looks odd! Calculate the PDF: \pause \[ h(x) = \tan(x), \quad h^{-1}(y) = \tan^{-1}(y). \] \pause \begin{eqnarray*} f_{Y}(y) &=& \frac{1}{\pi} \frac{d}{dy}\left(\tan^{-1} y \right) \quad -\infty < y < \infty, \pause \\ &=&\frac{1}{\pi(1+y^{2})} \quad -\infty