\documentclass[12pt]{article} \usepackage{../Math203} \usepackage{amsmath,amssymb} \usepackage{hyperref} \usepackage{color} \usepackage{enumitem} \pagestyle{plain} \renewcommand{\labelenumi}{\arabic{enumi}.} \renewcommand{\labelenumii}{\arabic{enumi}.\arabic{enumii}} \renewcommand{\labelenumiii}{\arabic{enumi}.\arabic{enumii}.\arabic{enumiii}} \setlength{\topmargin}{-0.75 in} \setlength{\textheight}{9.5 in} \parindent0in \newenvironment{packed_enum}{ \begin{enumerate} \setlength{\itemsep}{12pt} \setlength{\parskip}{0pt} \setlength{\parsep}{0pt} \setlength{\topsep}{0pt} }{\end{enumerate}} \renewcommand{\labelitemii}{$\blacktriangleright$} \begin{document} \begin{center} \textsclarge{Math 203: Principles of Statistics 1} \medskip \textsc{Fall 2024} \end{center} \begin{center} \textsc{Histograms} \end{center} <>= library(knitr) # global chunk options opts_chunk$set(cache=FALSE, autodep=TRUE) options(scipen=5) options(repos=c(CRAN="https://cloud.r-project.org/")) options(width=75) @ <>= myrnd <- function(x,d) {formatC(x, format="f", digits=d)} library(shape) library(diagram) @ The following data represent the birthweights (in grammes) of $n=189$ newborns. <>= library(MASS) x<-birthwt$bwt x @ We can plot a simple histogram: the default plot in \texttt{R} looks like this. <>= hist(x) freqs<-hist(x,plot=FALSE)$counts brks<-hist(x,plot=FALSE)$breaks @ We can see from the plot that there are $k=\Sexpr{length(brks)-1}$ intervals, and that the interval end points are \[ \Sexpr{brks} \] but we cannot easily see the individual frequencies for each category; they are \[ \Sexpr{freqs}. \] The default plot in \texttt{R} does not use the categorization rule as described in lectures, that is, the default rule in \texttt{R} puts a data value $x$ into an interval with endpoints $a$ and $b$ if \[ a < x \leq b \] whereas the method outlined in lectures uses the rule \[ a \leq x < b. \] Here, however, if we use the rule from lectures, we get <>= hist(x,right=FALSE) freqs<-hist(x,right=FALSE,plot=FALSE)$counts brks<-hist(x,right=FALSE,plot=FALSE)$breaks @ We have the individual frequencies for each category \[ \Sexpr{freqs} \] and there is no difference, as no individual data point matches any of the interval end points. \medskip We can add plotting options to make the plot more interpretable. In the next plot the following aspects have been altered: \begin{itemize} \item the number of intervals is increased to 18 covering the same range, \item the $x$-axis is re-labelled to be more informative, \item an informative title has been added. \end{itemize} <>= hist(x,breaks=seq(500,5000,by=250),right=FALSE,main='Birthweight data for n=189 newborns', xlab='Birthweight (g)') box() freqs<-hist(x,breaks=seq(500,5000,by=250),right=FALSE,plot=FALSE)$counts brks<-hist(x,breaks=seq(500,5000,by=250),right=FALSE,plot=FALSE)$breaks @ In the next plot, the same data are used in a \textbf{density} plot: <>= hist(x,breaks=seq(500,5000,by=250),right=FALSE,main='Birthweight data for n=189 newborns', xlab='Birthweight (g)',freq=FALSE) box() @ Note the change in the vertical scale: the frequency plot uses the frequencies \[ \Sexpr{freqs} \] to determine the height of each rectangle, whereas the density version takes a relative frequency, and divides it by the interval width (the `density' is measured per unit of the $x$-axis). In this plot the interval width is 250, so the density value for the eighth interval is \[ \frac{\textrm{Frequency in eighth interval}}{\textrm{Total sample size} \times \textrm{Interval width}} = \frac{\Sexpr{freqs[8]}}{\Sexpr{sum(freqs)} \times 250} = \Sexpr{myrnd(freqs[8]/(189*250),4)} \] In this case, there is no difference in the appearance of the plot, simply a change in the vertical axis values. This is because all the interval widths are the same. \bigskip We can also plot a \textbf{relative frequency} plot: <>= library(lattice) histogram(x,breaks=seq(500,5000,by=250),right=FALSE,main='Birthweight data for n=189 newborns', xlab='Birthweight (g)',col='gray') freqs<-hist(x,breaks=seq(500,5000,by=250),right=FALSE,plot=FALSE)$counts brks<-hist(x,breaks=seq(500,5000,by=250),right=FALSE,plot=FALSE)$breaks @ Here the relative frequencies are the frequencies per interval divided by the sample size. The relative frequencies are interval frequencies listed above divided by $n=189$ (here rounded to 4 decimal places): <>= round(freqs/189,4) @ In the plot, these relative frequencies are plotted as percentages (i.e. multiplied by 100). <>= round(100*freqs/189,4) @ \pagebreak It is possible to use unequal-sized intervals, but this makes interpretation of the histogram slightly more complicated. Suppose we want to used the intervals \[ [500,2000), [2000,2250), [2250,2500), [2500,2750), [2750,3000), [3250,3500), [3500,4000), [4000,5000] \] that is with the first interval width 1500g, the eighth interval width being 500g, the last interval width being 1000g, and the other interval widths being 250g. The corresponding histogram is <>= br<-c(500,seq(2000,3500,by=250),4000,5000) hist(x,breaks=br,right=FALSE,main='Birthweight data for n=189 newborns', xlab='Birthweight (g)',freq=TRUE) box() freqs<-hist(x,breaks=br,right=FALSE,plot=FALSE)$counts brks<-hist(x,breaks=br,right=FALSE,plot=FALSE)$breaks @ There are $k=\Sexpr{length(brks)-1}$ intervals and the individual frequencies for each of the intervals are \[ \Sexpr{freqs}. \] Due to the fact that the intervals are unequal sizes, it is preferable to use a density plot: here the density value for each interval is computed using the interval width for that interval. For example, for the first interval, we have \[ \frac{\textrm{Frequency in first interval}}{\textrm{Total sample size} \times \textrm{Interval width for first interval}} = \frac{\Sexpr{freqs[1]}}{\Sexpr{sum(freqs)} \times \Sexpr{brks[2]-brks[1]}} = \Sexpr{myrnd(freqs[1]/(189*(brks[2]-brks[1])),5)} \] whereas for the fifth interval, we have \[ \frac{\textrm{Frequency in fifth interval}}{\textrm{Total sample size} \times \textrm{Interval width for fifth interval}} = \frac{\Sexpr{freqs[5]}}{\Sexpr{sum(freqs)} \times \Sexpr{brks[6]-brks[5]}} = \Sexpr{myrnd(freqs[5]/(189*(brks[6]-brks[5])),5)} \] <>= br<-c(500,seq(2000,3500,by=250),4000,5000) hist(x,breaks=br,right=FALSE,main='Birthweight data for n=189 newborns', xlab='Birthweight (g)') box() @ \bigskip \begin{center} \textsc{Summary} \end{center} Key elements of the histogram: \begin{enumerate}[label=\arabic*.] \item \textbf{Number of intervals:} how many rectangles are used ? \item \textbf{Interval boundaries:} \begin{itemize} \item do the rectangles have equal widths or unequal widths ? \item which of the rules \[ a \leq x < b \qquad \textrm{or} \qquad a < x \leq b \] is being used to categorize data points lying exactly on interval boundaries ? \end{itemize} \item \textbf{Horizontal scale and label:} what is the variable of interest and its measurement unit ? \item \textbf{Vertical scale:} what type of vertical scale is being used \begin{itemize} \item frequency ? \item relative frequency ? \item density ? \end{itemize} \item \textbf{Title.} \end{enumerate} \end{document}