% Stats intro for Applied Stat I % Notes and comments at the end % \documentclass[serif]{beamer} % Serif for Computer Modern math font. \documentclass[serif, handout]{beamer} % Handout mode to ignore pause statements \hypersetup{colorlinks,linkcolor=,urlcolor=red} \usefonttheme{serif} % Looks like Computer Modern for non-math text -- nice! \setbeamertemplate{navigation symbols}{} % Supress navigation symbols % \usetheme{Berlin} % Displays sections on top \usepackage[english]{babel} % \definecolor{links}{HTML}{2A1B81} % \definecolor{links}{red} \setbeamertemplate{footline}[frame number] \mode % \mode{\setbeamercolor{background canvas}{bg=black!5}} \title{A Frequentist Introduction\footnote{See last slide for copyright information.}} \subtitle{STA442/2101 Fall 2018} \date{} % To suppress date \begin{document} \begin{frame} \titlepage \end{frame} \begin{frame} \frametitle{Background Reading} \framesubtitle{Optional} {\Large Chapter 1 of Davison's \emph{Statistical models}: Data, and probability models for data. } \end{frame} \begin{frame} \frametitle{Goal of statistical analysis} {\Large The goal of statistical analysis is to draw reasonable conclusions from noisy numerical data. } \end{frame} \begin{frame}{Steps in the process of statistical analysis} {One approach} \begin{itemize} \item Consider a fairly realistic example or problem. \pause \item Decide on a statistical model. \pause \item Perhaps decide sample size. \pause \item Acquire data. \pause \item Examine and clean the data; generate displays and descriptive statistics. \pause \item Estimate model parameters, for example by maximum likelihood. \pause \item Carry out tests, compute confidence intervals, or both. \pause \item Perhaps re-consider the model and go back to estimation. \pause \item Based on the results of estimation and inference, draw conclusions about the example or problem. \end{itemize} \end{frame} \begin{frame} \frametitle{What is a statistical model?} \framesubtitle{You should always be able to state the model.} \pause A \emph{statistical model} is a set of assertions that partly specify the probability distribution of the observable data. The specification may be direct or indirect. \pause \begin{itemize} \item Let $X_1, \ldots, X_n$ be a random sample from a normal distribution with expected value $\mu$ and variance $\sigma^2$. \pause \linebreak The parameters $\mu$ and $\sigma^2$ are unknown. \pause \item For $i=1, \ldots, n$, let $y_i = \beta_0 + \beta_1 x_{i,1} + \cdots + \beta_{p-1} x_{i,p-1} + \epsilon_i$, where \begin{itemize} \item[] $\beta_0, \ldots, \beta_{p-1}$ are unknown constants. \item[] $x_{i,j}$ are known constants. \item[] $\epsilon_1, \ldots, \epsilon_n$ are independent $N(0,\sigma^2)$ random variables. \item[] $\sigma^2$ is an unknown constant. \item[] $y_1, \ldots, y_n$ are observable random variables. \pause \end{itemize} The parameters $\beta_0, \ldots, \beta_{p-1}, \sigma^2$ are unknown. \end{itemize} \end{frame} \begin{frame} \frametitle{Model and Truth} \framesubtitle{Is a statistical model the same thing as the truth?} \pause \begin{quote} ``Essentially all models are wrong, but some are useful." (Box and Draper, 1987, p. 424) \end{quote} % Box, G. E. P. and Draper, N. R. (1987). Empirical Model-Building and Response Surfaces. New York: Wiley. % \vspace{5mm} \pause % It helps if the unknown parameters represent what you wish you knew about the data. % It also helps if the statistical model is a halfway believable substantive model for how the data might have been produced, and if the distribution implied by the model is similar to the empirical distribution of the data. \end{frame} \begin{frame} \frametitle{Parameter Space} The \emph{parameter space} is the set of values that can be taken on by the parameter. \pause \begin{itemize} \item Let $X_1, \ldots, X_n$ be a random sample from a normal distribution with expected value $\mu$ and variance $\sigma^2$. \pause The parameter space is $\{(\mu,\sigma^2): -\infty < \mu < \infty, \sigma^2 > 0\}$. \pause \item For $i=1, \ldots, n$, let $y_i = \beta_0 + \beta_1 x_{i,1} + \cdots + \beta_{p-1} x_{i,p-1} + \epsilon_i$, where \begin{itemize} \item[] $\beta_0, \ldots, \beta_{p-1}$ are unknown constants. \item[] $x_{i,j}$ are known constants. \item[] $\epsilon_1, \ldots, \epsilon_n$ are independent $N(0,\sigma^2)$ random variables. \item[] $\sigma^2$ is an unknown constant. \item[] $y_1, \ldots, y_n$ are observable random variables. \end{itemize} \pause The parameter space is $\{(\beta_0, \ldots, \beta_{p-1}, \sigma^2): -\infty < \beta_j < \infty, \sigma^2 > 0\}$. \end{itemize} \end{frame} \begin{frame}{Coffee taste test} A fast food chain is considering a change in the blend of coffee beans they use to make their coffee. To determine whether their customers prefer the new blend, the company plans to select a random sample of $n=100$ coffee-drinking customers and ask them to taste coffee made with the new blend and with the old blend, in cups marked ``$A$" and ``$B$." Half the time the new blend will be in cup $A$, and half the time it will be in cup $B$. Management wants to know if there is a difference in preference for the two blends. \end{frame} \begin{frame}{Statistical model} Letting $\theta$ denote the probability that a consumer will choose the new blend, treat the data $Y_1, \ldots, Y_n$ as a random sample from a Bernoulli distribution. That is, independently for $i=1, \ldots, n$, \begin{displaymath} P(y_i|\theta) = \theta^{y_i} (1-\theta)^{1-y_i} \end{displaymath} for $y_i=0$ or $y_i=1$, and zero otherwise. % The conditional probability notation is not in the book (I believe). \vspace{5mm} \pause \begin{itemize} \item Parameter space is the interval from zero to one. \pause \item $\theta$ could be estimated by maximum likelihood. \pause \item Large-sample tests and confidence intervals are available. \end{itemize} \pause Note that $Y = \sum_{i=1}^n Y_i$ is the number of consumers who choose the new blend. Because $Y \sim B(n,\theta)$, the whole experiment could also be treated as a single observation from a Binomial. \end{frame} \begin{frame}{Find the MLE of $\theta$}{Show your work} \pause Denoting the likelihood by $L(\theta)$ and the log likelihood by $\ell(\theta) = \log L(\theta)$, maximize the log likelihood. \pause \begin{eqnarray*} \frac{\partial\ell}{\partial\theta} \pause & = & \frac{\partial}{\partial\theta} \log\left(\prod_{i=1}^n P(y_i|\theta) \right) \\ \pause & = & \frac{\partial}{\partial\theta} \log\left(\prod_{i=1}^n \theta^{y_i} (1-\theta)^{1-y_i} \right) \\ \pause & = & \frac{\partial}{\partial\theta} \log\left(\theta^{\sum_{i=1}^n y_i} (1-\theta)^{n-\sum_{i=1}^n y_i}\right) \\ \pause & = & \frac{\partial}{\partial\theta}\left((\sum_{i=1}^n y_i)\log\theta + (n-\sum_{i=1}^n y_i)\log (1-\theta) \right) \\ \pause & = & \frac{\sum_{i=1}^n y_i}{\theta} - \frac{n-\sum_{i=1}^n y_i}{1-\theta} \end{eqnarray*} \end{frame} \begin{frame}{Setting the derivative to zero and solving} \pause \begin{itemize} \item $\theta = \frac{\sum_{i=1}^n y_i}{n} = \overline{y}$ \pause \item Second derivative test: $ \frac{\partial^2\log \ell}{\partial\theta^2} = -n\left(\frac{1-\overline{y}}{(1-\theta)^2} + \frac{\overline{y}}{\theta^2} \right) < 0$ \pause \item Concave down, maximum, and the MLE is the sample proportion: \pause $\widehat{\theta} = \overline{y} = p$ \end{itemize} \end{frame} \begin{frame}[fragile] % Note use of fragile to make verbatim work. \frametitle{Numerical estimate} Suppose 60 of the 100 consumers prefer the new blend. Give a point estimate the parameter $\theta$. Your answer is a number. \vspace{10mm} \pause \begin{verbatim} > p = 60/100; p [1] 0.6 \end{verbatim} \end{frame} \begin{frame} \frametitle{Tests of statistical hypotheses} %\framesubtitle{} \pause \begin{itemize} \item Model: $Y \sim F_\theta $ \pause \item $Y$ is the data vector, and $\mathcal{Y}$ is the sample space: $Y \in \mathcal{Y}$ \pause \item $\theta$ is the parameter, and $\Theta$ is the parameter space: $\theta \in \Theta$ \pause \item Null hypothesis is $H_0: \theta \in \Theta_0 \mbox{ v.s. } H_A: \theta \in \Theta \cap \Theta_0^c$. \pause \item Meaning of the \emph{null} hypothesis is that \emph{nothing} interesting is happening. \pause \item $\mathcal{C} \subset \mathcal{Y}$ is the \emph{critical region}. Reject $H_0$ in favour of $H_A$ when $Y \in \mathcal{C}$. \pause \item Significance level $\alpha$ (\emph{size} of the test) is the maximum probability of rejecting $H_0$ when $H_0$ is true. \pause Conventionally, $\alpha=0.05$. \pause \item $p$-value is the smallest value of $\alpha$ for which $H_0$ can be rejected. \pause \item Small $p$-values are interpreted as providing stronger evidence against the null hypothesis. \end{itemize} \end{frame} \begin{frame} \frametitle{Type I and Type II error} \framesubtitle{A Neyman-Pearson idea rather than Fisher} \pause \begin{itemize} \item Type I error is to reject $H_0$ when $H_0$ is true. \pause \item Type II error is to \emph{not} reject $H_0$ when $H_0$ is false. \pause % \item Can't minimize the probability of both types of error at once. % \item So hold the maximum probability of a Type I error (the significance level $\alpha$) to a ``small" value like $\alpha=0.05$, and then seek tests that minimize the probability of Type II error. \item $1-Pr\{$Type II Error$\}$ is called \emph{power}. \pause \item If two tests have the same maximum Type I error probability $\alpha$, the one with higher power is better. \pause \item Power may also be used to select sample size. \end{itemize} \end{frame} \begin{frame} \frametitle{Carry out a test to determine which brand of coffee is preferred} \framesubtitle{Recall the model is $Y_1, \ldots, Y_n \stackrel{i.i.d.}{\sim} B(1,\theta)$} \pause Start by stating the null hypothesis. \pause \begin{itemize} \item $H_0: \theta=0.50$ \item $H_1: \theta \neq 0.50$ \pause \item Could you make a case for a one-sided test? \pause \item $\alpha=0.05$ as usual. \pause \item Central Limit Theorem says $\widehat{\theta}=\overline{Y}$ is approximately normal with mean $\theta$ and variance $\frac{\theta(1-\theta)}{n}$. \end{itemize} \end{frame} \begin{frame}[fragile] \frametitle{Several valid test statistics for $H_0: \theta=\theta_0$ are available} \framesubtitle{Recall that approximately, $\overline{Y} \sim N(\theta,\frac{\theta(1-\theta)}{n})$} \pause Two of them are % Which one do you like more? Why? \begin{displaymath} Z_1 = \frac{\sqrt{n}(\overline{Y}-\theta_0)}{\sqrt{\theta_0(1-\theta_0)}} \end{displaymath} and \pause \begin{displaymath} Z_2 = \frac{\sqrt{n}(\overline{Y}-\theta_0)}{\sqrt{\overline{Y}(1-\overline{Y})}} \end{displaymath} \vspace{10mm} \pause What is the critical value? Your answer is a number. \pause \begin{verbatim} > alpha = 0.05 > qnorm(1-alpha/2) [1] 1.959964 \end{verbatim} \end{frame} \begin{frame}[fragile] \frametitle{Calculate the test statistic and the $p$-value for each test} \framesubtitle{Suppose 60 out of 100 preferred the new blend} \pause $ Z_1 = \frac{\sqrt{n}(\overline{Y}-\theta_0)}{\sqrt{\theta_0(1-\theta_0)}}$ \pause \begin{verbatim} > theta0 = .5; ybar = .6; n = 100 > Z1 = sqrt(n)*(ybar-theta0)/sqrt(theta0*(1-theta0)); Z1 [1] 2 > pval1 = 2 * (1-pnorm(Z1)); pval1 [1] 0.04550026 \end{verbatim} \pause $Z_2 = \frac{\sqrt{n}(\overline{Y}-\theta_0)}{\sqrt{\overline{Y}(1-\overline{Y})}}$ \pause \begin{verbatim} > Z2 = sqrt(n)*(ybar-theta0)/sqrt(ybar*(1-ybar)); Z2 [1] 2.041241 > pval2 = 2 * (1-pnorm(Z2)); pval2 [1] 0.04122683 \end{verbatim} \end{frame} \begin{frame} \frametitle{Conclusions} %\framesubtitle{In symbols and words: Words are more important} \begin{itemize} \item Do you reject $H_0$? \pause \emph{Yes, just barely.} \pause \item Isn't the $\alpha=0.05$ significance level pretty arbitrary? \pause \linebreak \emph{Yes, but if people insist on a Yes or No answer, this is what you give them.} \pause \item What do you conclude, in symbols? \pause $\theta \neq 0.50$. \emph{Specifically,} $\theta > 0.50$. \pause \item What do you conclude, in plain language? Your answer is a statement about coffee. \pause \emph{More consumers prefer the new blend of coffee beans.} \pause \item Can you really draw directional conclusions when all you did was reject a non-directional null hypothesis? \pause \emph{Yes.} \end{itemize} \end{frame} \begin{frame} \frametitle{A technical issue} %\framesubtitle{} { \small \begin{itemize} \item In this class we will mostly avoid one-tailed tests. \pause \item Why? Ask what would happen if the results were strong and in the opposite direction to what was predicted (dental example). \pause \item But when $H_0$ is rejected, we still draw directional conclusions. \pause \item For example, if $x$ is income and $y$ is credit card debt, we test $H_0: \beta_1=0$ with a two-sided $t$-test. \pause \item Say $p = 0.0021$ and $\widehat{\beta}_1 = 1.27$. \pause We say ``Consumers with higher incomes tend to have more credit card debt." \pause \item Is this justified? We'd better hope so, or all we can say is ``There is a connection between income and average credit card debt." \pause \item Then they ask: ``What's the connection? Do people with lower income have more debt?" \pause \item And you have to say ``Sorry, I don't know." \pause \item It's a good way to get fired, or at least look silly. \end{itemize} } % End size \end{frame} \begin{frame} \frametitle{The technical resolution} %\framesubtitle{} Decompose the two-sided test into a set of two one-sided tests with significance level $\alpha/2$, equivalent to the two-sided test. \end{frame} \begin{frame} \frametitle{Two-sided test} %\framesubtitle{} \begin{center} {\Large $H_0: \theta=\frac{1}{2}$ versus $H_1: \theta \neq \frac{1}{2}$, $\alpha=0.05$ } \vspace{10mm} \includegraphics[width=4.5in]{bothtails} \end{center} \end{frame} \begin{frame} \frametitle{Left-sided test} %\framesubtitle{} \begin{center} {\Large $H_0: \theta\geq \frac{1}{2}$ versus $H_1: \theta < \frac{1}{2}$, $\alpha=0.05$ } \vspace{10mm} \includegraphics[width=4.5in]{lefttail} \end{center} \end{frame} \begin{frame} \frametitle{Right-sided test} %\framesubtitle{} \begin{center} {\Large $H_0: \theta\leq \frac{1}{2}$ versus $H_1: \theta > \frac{1}{2}$, $\alpha=0.05$ } \vspace{10mm} \includegraphics[width=4.5in]{righttail} \end{center} \end{frame} \begin{frame} \frametitle{Decomposing the 2-sided test into two 1-sided tests} \pause %\framesubtitle{} \begin{center} \begin{tabular}{lc} \raisebox{0.25in}{\small $H_0: \theta=\frac{1}{2}$ vs. $H_1: \theta \neq \frac{1}{2}$, $\alpha=0.05$} & \includegraphics[width=2in]{bothtails} \\ \raisebox{0.25in}{\small $H_0: \theta\geq \frac{1}{2}$ vs. $H_1: \theta < \frac{1}{2}$, $\alpha=0.05$} & \includegraphics[width=2in]{lefttail} \\ \raisebox{0.25in}{\small $H_0: \theta\leq \frac{1}{2}$ versus $H_1: \theta > \frac{1}{2}$, $\alpha=0.05$} & \includegraphics[width=2in]{righttail} \\ \end{tabular} \end{center} \pause \begin{itemize} \item Clearly, the 2-sided test rejects $H_0$ if and only if exactly \emph{one} of the 1-sided tests reject $H_0$. \pause \item Carry out \emph{both} of the one-sided tests. \pause \item Draw a directional conclusion if $H_0$ is rejected. \end{itemize} \end{frame} \begin{frame} \frametitle{Summary of the technical resolution} \pause %\framesubtitle{} \begin{itemize} \item Decompose the two-sided test into a set of two one-sided tests with significance level $\alpha/2$, equivalent to the two-sided test. \pause \item In practice, just look at the sign of the regression coefficient, or compare the sample means. \pause \item Under the surface you are decomposing the two-sided test, but you never mention it. \end{itemize} \end{frame} \begin{frame} \frametitle{Plain language} \pause %\framesubtitle{} \begin{itemize} \item It is very important to state directional conclusions, and state them clearly in terms of the subject matter. \textbf{Say what happened!} If you are asked state the conclusion in plain language, your answer \emph{must} be free of statistical mumbo-jumbo. \pause \item \emph{Marking rule}: If the question asks for plain language and you draw a non-directional conclusion when a directional conclusion is possible, you get half marks at most. \end{itemize} \end{frame} \begin{frame} \frametitle{What about negative conclusions?} \framesubtitle{What would you say if $Z=1.84$?} \pause Here are two possibilities, in plain language. \pause \begin{itemize} \item ``This study does not provide clear evidence that consumers prefer one blend of coffee beans over the other." \pause \item ``The results are consistent with no difference in preference for the two coffee bean blends." \end{itemize} \vspace{5mm} \pause In this course, we will not just casually accept the null hypothesis. \pause We will \emph{not} say that there was no difference in preference. \pause \vspace{3mm} We are taking the side of Fisher over Neyman and Pearson in an old and very nasty philosophic dispute. \end{frame} \begin{frame} \frametitle{Confidence intervals} \framesubtitle{Usually for individual parameters} \pause \begin{itemize} \item Point estimates may give a false sense of precision. \pause \item We should provide a margin of probable error as well. \end{itemize} \end{frame} \begin{frame} \frametitle{Confidence Intervals} \framesubtitle{Taste test example} \pause Approximately for large $n$, \pause \begin{eqnarray*} 1-\alpha & = & Pr\{ -z_{\alpha/2} < Z < z_{\alpha/2} \} \\ \pause & \approx & Pr\left\{ -z_{\alpha/2} < \frac{\sqrt{n}(\overline{Y}-\theta)}{\sqrt{\overline{Y}(1-\overline{Y})}} < z_{\alpha/2} \right\} \\ \pause & = & Pr\left\{ \overline{Y} - z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{n}} < \theta < \overline{Y} + z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{n}} \right\} \end{eqnarray*} \pause \begin{itemize} \item Could express this as $\overline{Y} \pm z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{n}}$. \pause \item $z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{n}}$ is sometimes called the \emph{margin of error}. \pause \item If $\alpha=0.05$, it's the 95\% margin of error. \end{itemize} \end{frame} \begin{frame} \frametitle{Give a 95\% confidence interval for the taste test data.} \framesubtitle{The answer is a pair of numbers. Show some work.} \pause \begin{eqnarray*} & & \left(\overline{y} - z_{\alpha/2}\sqrt{\frac{\overline{y}(1-\overline{y})}{n}} ~~,~~ \overline{y} + z_{\alpha/2}\sqrt{\frac{\overline{y}(1-\overline{y})}{n}} \right) \\ & & \\ &=& \left(0.60 - 1.96\sqrt{\frac{0.6\times 0.4}{100}} ~~,~~ 0.60 + 1.96\sqrt{\frac{0.6\times 0.4}{100}} \right) \\ & & \\ &=& (0.504,0.696) \end{eqnarray*} \pause In a report, you could say \pause \begin{itemize} \item The estimated proportion preferring the new coffee bean blend is $0.60 \pm 0.096$, or \pause \item ``Sixty percent of consumers preferred the new blend. These results are expected to be accurate within 10 percentage points, 19 times out of 20." \end{itemize} \end{frame} \begin{frame} \frametitle{Meaning of the confidence interval} \begin{itemize} \item We calculated a 95\% confidence interval of $(0.504,0.696)$ for $\theta$. \item Does this mean $Pr\{ 0.504 < \theta < 0.696 \}=0.95$? \pause \item No! The quantities $0.504$, $0.696$ and $\theta$ are all constants, so $Pr\{ 0.504 < \theta < 0.696 \}$ is either zero or one. \pause \item The endpoints of the confidence interval are random variables, and the numbers $0.504$ and $0.696$ are \emph{realizations} of those random variables, arising from a particular random sample. \pause \item Meaning of the probability statement: If we were to calculate an interval in this manner for a large number of random samples, the interval would contain the true parameter around $95\%$ of the time. \pause \item The confidence interval is a guess, and the guess is either right or wrong. But the guess is the constructed by a method that is right 95\% of the time. % Take it or leave it. \end{itemize} \end{frame} \begin{frame} \frametitle{More on confidence intervals} %\framesubtitle{Usually, confidence \emph{intervals} for single parameters} \begin{itemize} \item Can have confidence \emph{regions} for the entire parameter vector or multi-dimensional functions of the parameter vector. \item Confidence regions correspond to tests. \end{itemize} \end{frame} \begin{frame} \frametitle{Confidence intervals (regions) correspond to tests} \framesubtitle{Recall $Z_1 = \frac{\sqrt{n}(\overline{Y}-\theta_0)} {\sqrt{\theta_0(1-\theta_0)}}$ and $Z_2 = \frac{\sqrt{n} (\overline{Y}-\theta_0)}{\sqrt{\overline{Y}(1-\overline{Y})}}$.} \pause $H_0$ is \emph{not} rejected if and only if \begin{displaymath} -z_{\alpha/2} < Z_2 < z_{\alpha/2} \end{displaymath} \pause if and only if \begin{displaymath} \overline{Y} - z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{n}} < \theta_0 < \overline{Y} + z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{n}} \end{displaymath} \pause % From the derivation of the confidence interval, \begin{itemize} \item So the confidence interval consists of those parameter values $\theta_0$ for which $H_0: \theta=\theta_0$ is \emph{not} rejected. \pause \item That is, the null hypothesis is rejected at significance level $\alpha$ if and only if the value given by the null hypothesis is outside the $(1-\alpha)\times 100\%$ confidence interval. % \item There is a confidence interval corresponding to $Z_1$ too. % \item In general, any test can be inverted to obtain a confidence region. \end{itemize} \end{frame} \begin{frame} \frametitle{Selecting sample size} \begin{itemize} \item Where did that $n=100$ come from? \pause \item Probably off the top of someone's head. \pause \item We can (and should) be more systematic. \pause \item Sample size can be selected \begin{itemize} \item To achieve a desired margin of error \item To achieve a desired statistical power \item In other reasonable ways \end{itemize} \end{itemize} \end{frame} \begin{frame} \frametitle{Statistical Power} \pause The power of a test is the probability of rejecting $H_0$ when $H_0$ is false. \pause \begin{itemize} \item More power is good. \pause \item Power is not just one number. It is a \emph{function} of the parameter(s). \pause \item Usually, \begin{itemize} \item For any $n$, the more incorrect $H_0$ is, the greater the power. \pause \item For any parameter value satisfying the alternative hypothesis, the larger $n$ is, the greater the power. \end{itemize} \end{itemize} \end{frame} \begin{frame} \frametitle{Statistical power analysis} \framesubtitle{To select sample size} \pause \begin{itemize} \item Pick an effect you'd like to be able to detect -- a parameter value such that $H_0$ is false. It should be just over the boundary of interesting and meaningful. \pause \item Pick a desired power, a probability with which you'd like to be able to detect the effect by rejecting the null hypothesis. \pause \item Start with a fairly small $n$ and calculate the power. Increase the sample size until the desired power is reached. \pause \end{itemize} \vspace{5mm} There are two main issues. \vspace{5mm} \begin{itemize} \item What is an ``interesting" or ``meaningful" parameter value? \pause \item How do you calculate the probability of rejecting $H_0$? \end{itemize} \end{frame} \begin{frame} \frametitle{Calculating power for the test of a single proportion} \framesubtitle{True parameter value is $\theta$} \pause %\begin{displaymath} % Z_2 = \frac{\sqrt{n}(\overline{Y}-\theta_0)}{\sqrt{\overline{Y}(1-\overline{Y})}} %\end{displaymath} {\tiny \begin{eqnarray*} \mbox{Power} &=& 1-Pr\{-z_{\alpha/2} < Z_2 < z_{\alpha/2} \} \\ \pause &=& 1-Pr\left\{-z_{\alpha/2} < \frac{\sqrt{n}(\overline{Y}-\theta_0)}{\sqrt{\overline{Y}(1-\overline{Y})}} < z_{\alpha/2} \right\} \\ \pause && \\ &=& \ldots \\ \pause && \\ &=& 1-Pr\left\{\frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} - z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{\theta(1-\theta)}} ~~~<~~~ {\color{red}\frac{\sqrt{n}(\overline{Y}-\theta)}{\sqrt{\theta(1-\theta)}} } \right. \\ & & \hspace{8mm} < \left. \frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} + z_{\alpha/2}\sqrt{\frac{\overline{Y}(1-\overline{Y})}{\theta(1-\theta)}} \right\}\\ \pause &\approx& 1-Pr\left\{\frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} - z_{\alpha/2} < { \color{red} Z } < \frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} + z_{\alpha/2} \right\} \\ \pause &=& 1 - \Phi\left( \frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} + z_{\alpha/2} \right) + \Phi\left( \frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} - z_{\alpha/2} \right), \end{eqnarray*} } \pause where $\Phi(\cdot)$ is the cumulative distribution function of the standard normal. \end{frame} \begin{frame}[fragile] \frametitle{An R function to calculate approximate power} \framesubtitle{For the test of a single proportion} \pause {\small \begin{displaymath} \mbox{Power} = 1 - \Phi\left( \frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} + z_{\alpha/2} \right) + \Phi\left( \frac{\sqrt{n}(\theta_0-\theta)}{\sqrt{\theta(1-\theta)}} - z_{\alpha/2} \right) \end{displaymath} } \pause {\small \begin{verbatim} Z2power = function(theta,n,theta0=0.50,alpha=0.05) { effect = sqrt(n)*(theta0-theta)/sqrt(theta*(1-theta)) z = qnorm(1-alpha/2) Z2power = 1 - pnorm(effect+z) + pnorm(effect-z) Z2power } # End of function Z2power \end{verbatim} } \end{frame} \begin{frame}[fragile] \frametitle{Some numerical examples} \framesubtitle{\texttt{Z2power = function(theta,n,theta0=0.50,alpha=0.05)}} \pause \begin{verbatim} > Z2power(0.50,100) # Should be alpha = 0.05 [1] 0.05 > > Z2power(0.55,100) [1] 0.1713209 > Z2power(0.60,100) [1] 0.5324209 > Z2power(0.65,100) [1] 0.8819698 > Z2power(0.40,100) [1] 0.5324209 > Z2power(0.55,500) [1] 0.613098 > Z2power(0.55,1000) [1] 0.8884346 \end{verbatim} \end{frame} \begin{frame}[fragile] \frametitle{Find smallest sample size needed to detect $\theta=0.60$ as different from $\theta_0=0.50$ with probability at least 0.80} \pause \begin{verbatim} > samplesize = 1 > power=Z2power(theta=0.60,n=samplesize); power [1] 0.05478667 > while(power < 0.80) + { + samplesize = samplesize+1 + power = Z2power(theta=0.60,n=samplesize) + } > samplesize [1] 189 > power [1] 0.8013024 \end{verbatim} \end{frame} \begin{frame} \frametitle{What is required of the scientist} \framesubtitle{Who wants to select sample size by power analysis} \pause The scientist must specify \pause \begin{itemize} \item Parameter values that he or she wants to be able to detect as different from $H_0$ value. \pause \item Desired power (probability of detection) \end{itemize} \vspace{20mm} \pause It's not always easy for a scientist to think in terms of the parameters of a statistical model. \end{frame} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{frame} \frametitle{Copyright Information} This slide show was prepared by \href{http://www.utstat.toronto.edu/~brunner}{Jerry Brunner}, Department of Statistics, University of Toronto. It is licensed under a \href{http://creativecommons.org/licenses/by-sa/3.0/deed.en_US} {Creative Commons Attribution - ShareAlike 3.0 Unported License}. Use any part of it as you like and share the result freely. The \LaTeX~source code is available from the course website: \href{http://www.utstat.toronto.edu/~brunner/oldclass/appliedf18} {\footnotesize \texttt{http://www.utstat.toronto.edu/$^\sim$brunner/oldclass/appliedf18}} \end{frame} \end{document} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{frame} \frametitle{} %\framesubtitle{} \begin{itemize} \item \item \item \end{itemize} \end{frame} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Notes and comments In 2013, Cut out non-central chisq part, added nice what's a model, parameter space Non-central chisq is in 2012 version. In 2014, based on a one-vs. 2-sample t-test rather than binomial. Version 2014a wandered back into the taste test example. In 2016, lifted material from 2012 to beef it up and prepare for power of the general linear test. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Making pictures for that decomposition of the 2-tailed test rm(list=ls()) z = seq(from=-3.5,to=3.5,by=0.01) Density = dnorm(z) # Produces bothtails.pdf plot(z,Density,type='l',xlab="",ylab="",axes="F") # Draw base line x1 = c(-3.5,3.5); y1 = c(0,0); lines(x1,y1) # Draw vertical lines ht = dnorm(1.96) x2 = c(-1.96,-1.95); y2 = c(0,ht); lines(x2,y2) x3 = c(1.96,1.95); y3 = c(0,ht); lines(x3,y3) # Place 0.025 text text(-2.15,.015,"0.025"); text(2.15,.015,"0.025") # Produces lefttail.pdf plot(z,Density,type='l',xlab="",ylab="",axes="F") # Draw base line x1 = c(-3.5,3.5); y1 = c(0,0); lines(x1,y1) # Draw vertical line ht = dnorm(1.96) x2 = c(-1.96,-1.95); y2 = c(0,ht); lines(x2,y2) # x3 = c(1.96,1.95); y3 = c(0,ht); lines(x3,y3) # Place 0.025 text text(-2.15,.015,"0.025"); # text(2.15,.015,"0.025") # Produces righttail.pdf plot(z,Density,type='l',xlab="",ylab="",axes="F") # Draw base line x1 = c(-3.5,3.5); y1 = c(0,0); lines(x1,y1) # Draw vertical line ht = dnorm(1.96) # x2 = c(-1.96,-1.95); y2 = c(0,ht); lines(x2,y2) x3 = c(1.96,1.95); y3 = c(0,ht); lines(x3,y3) # Place 0.025 text # text(-2.15,.015,"0.025") text(2.15,.015,"0.025")