Lab 1

\[ \text{estimator} = \frac{1}{N}\sum_{i=1}^N Y_i \qqtext{ where } Y_1 \ldots Y_N \qqtext{ is our sample } \]

sam.plot + geom_hline(yintercept=estimator(sam))

• The horizontal line shows the mean of my sample \(Y_1 \ldots Y_n\).
• That is, the fraction of times I head ‘I prefer B’ when I ask someone’s preference.

• People often talk as if both \(Y_i\) and \(y_j\) are the same kind of thing. e.g. people’s preferences.
  • That’s a source of a lot of confusion because they aren’t.
  • To make this clear, it helps to be a bit concrete.
• Let’s suppose we’re doing a poll by phone. The way we find out someone’s preference is to call them.
  • \(y_j\) is a person’s preference. A real person. e.g. the \(j\)th one in the phone book.
  • \(Y_i\) is an utterance. It’s what we hear on our \(i\)th call.
• In this class, no matter what the context, I’ll say ‘\(i\)th call’ and ‘\(j\)th’ person.
  • It could be a quality control application in a brewery.
  • I’m opening and pouring a few beers to check that we haven’t overcarbonated.
  • My \(j\)th person is a can and my \(i\)th call is a glass.
• As you can imagine, that’s not what people do in real life. And that’s ok.
  • Once you’re used to this stuff, you know what they mean.
  • And if you don’t, you can ask. You can say, e.g., do you mean ‘little \(y\)’ or ‘big \(Y\)’?

Why \(Y_i\) and \(y_j\)?

• This notation helps remind us of what’s similar and what’s different.
  • We use the same letter to communicate that each (big) \(Y\) is somebody’s (little) \(y\).
  • We capitalize to communicate that \(Y_1\) isn’t \(y_1\).
  • In this class, I like to count them out with different letters too.
    • \(Y_i\) for \(i\) from \(1\) to \(n\). \(y_j\) from \(j\) from \(1\) to \(m\)
    • The other stuff is common practice. This counting stuff less so.
• This same-but-capitalized thing is helpful in writing but not so great speaking.
  • It gets confusing because \(Y\) and \(y\) are usually pronounced the same.
  • Often, it’s clear from context what you mean.
    • Especially if you say the subscript.
    • \(Y_i\) and \(y_j\) sound like ‘why eye’ and ‘why j’
  • But when it’s not, we say ‘big \(Y\)’ and ‘little \(y\)’

• It takes two arguments — a list and a function — and calls the function on every element of the list.
  • Then gives you back a list — or actually a vector, hence ’_vec’ — of the results.
  • Lists and vectors being different in R. You can give it either.
• vs. writing a for loop, this saves us the ‘filing’
  • It allocates the output vector and puts each function value into it.
  • To make it a one-liner, you can use an anonymous function.

one.to.ten.again = 0:9 |> map_vec(\(x) x+1)
one.to.ten.again

• This is a built-in function most modern lanuages. Python, Julia, Scheme, Haskell, ML, etc.
  • And it’s always called ‘map’. Except in R. It’s built in, but it’s called ‘sapply’.
  • I tend to use the purrr version just because the name is better.
  • People speak the code above this way: ‘map [the function] plus.one over the vector zero through nine’
• We’re using it for a reason that’s common but not the ‘usual reason’ people do.
  • We don’t actually care about getting passed every element of the list.
  • We just want to repeat something — something random — a bunch of times.
  • So we’ll pass it a list with the same thing repeated a bunch of times.
• We’ll write a function that …
  • takes as input a population
  • samples the population
  • returns a sample mean—an estimate of the population mean.
• We’ll run it R times on our population by mapping our function over a list containing R copies the population.

pops = rep(list(pop), R)

• You may have noticed that I said map_vec is function of two arguments.
  • A list. Above, that’s 0:9.
  • A function. Above, that’s plus.one.
• You’d think I’d call it like this.

map_vec(0:9, plus.one)

• Instead, I called it like this.

0:9 |> map_vec(plus.one)

• This is called ‘piping’. We say the list 0:9 is ‘piped in’.
  • a |> f(...) is the same as f(a,...) for any function f and list of arguments ....
• We use piping for the same reason we might say …
  • ‘take the list zero to nine and map plus.one over it’ …
  • vs. ‘map plus.one over the list zero to nine’.
• It’s especially natural when you want to do a sequence of things to a list.
  • e.g., if we want to take the list zero to nine, add one to each element, then square each result.

square = function(x) { x^2 }
0:9 |> map_vec(plus.one) |> map_vec(square) # in piped form
map_vec(map_vec(0:9, plus.one), square)     # in unpiped form

estimate = NaN                                       
sam.plot = pop.plot + geom_hline(yintercept=estimate)    
sam.plot          

• Color in the dots you talked to.
  • If talked to the same one multiple times, draw in another dot.
• Then draw a horizontal line indicating your estimate: the mean of your sample \(Y_1 \ldots Y_n\).
• I’ve made the variable estimate NaN. That stands for Not a Number. Change it to a number.
  • ggplot will just ignore it until you do, but it’ll draw it in once you have.
  • That’ll give you a sense of how precisely you’ve drawn your line.
• You can use the built-in whiteboard tool to do draw with your stylus, mouse, trackpad, etc.
  • You should be able to just draw with a stylus. If you want to scroll and change tabs instead, use a finger.
  • With a mouse or trackpad, double-click to open it the tool.
  • While in whiteboard mode, you can’t scroll, type, etc. Click the ‘x’ on the bottom panel to get out.

with.replacement = function(pop, sampling.rate=1/2) { 
  N = round(m * sampling.rate) 
  J = sample(1:m, N, replace=TRUE)
  sam = pop[J,]
  estimator(sam)
}

estimator.samples =  pops |> map_vec(with.replacement)

ggplot(data.frame(sample=estimator.samples, rep=(1:R)/R)) + 
  geom_bar(aes(x=sample, y=after_stat(prop)), alpha=.3) +
  geom_point(aes(x=sample, y=rep), position=position_jitter(w=.01,h=0), color='purple', alpha=.3) + 
  geom_vline(xintercept=estimate) 

Q. Why are there four bars?

• Some of you may have called the same person twice. That feels a little wasteful.
  • If you were willing to make three calls, you could’ve gained more information.
  • People tend not to like this kind of sampling in real life.
• But there is something to like. The math is easy.
• We say the \(Y\)s in our sample are independent and identically distributed.
  • What you hear on one call tells you nothing about what you’ll hear on another if you’ve seen the population.
  • No two calls are different. The probability of each possible response is the same on all calls.
• The independent thing can feel counterintuitive because they could be calls to the same person.
  • The second time you wouldn’t have to ask. You’d know their response.
  • Q. Why—informally—does this not make our \(Y\)s dependent?
  • Because it’s knowledge about a person, not a call.
  • If you’ve seen the population already, you wouldn’t have to ask anyone.

# editor-code-line-numbers: 2
without.replacement = function(pop, sampling.rate=1/2) {
  N = round(m * sampling.rate) 
  J = sample(1:m, N, replace=FALSE)
  sam = pop[J,]
  estimator(sam)
}

estimator.samples =  pops |> map_vec(without.replacement)

ggplot(data.frame(sample=estimator.samples, rep=(1:R)/R)) + 
  geom_bar(aes(x=sample, y=after_stat(prop)), alpha=.3) +
  geom_point(aes(x=sample, y=rep), position=position_jitter(w=.01,h=0), color='purple', alpha=.3) + 
  geom_vline(xintercept=estimate) 

In R, this a one line change to the sampling code.

• People like this because it doesn’t feel wasteful. You don’t call the same person twice.
• But it makes the math a bit harder. Calls aren’t independent.
• So people often do this, and do the math as if they’d sampled with replacement.
• The justification is that, with enough balls, you weren’t going to pick the same one twice anyway.
• What do you think of this? Do you think it makes sense?

# editor-code-line-numbers: 2-3
by.randomization = function(pop, sampling.rate=1/2) {
  flips = rbinom(m, 1, sampling.rate)
  sam = pop[flips==1,]
  estimator(sam)
}

estimator.samples =  pops |> map_vec(by.randomization)

ggplot(data.frame(sample=estimator.samples, rep=(1:R)/R)) + 
  geom_bar(aes(x=sample, y=after_stat(prop)), alpha=.3) +
  geom_point(aes(x=sample, y=rep), position=position_jitter(w=.01,h=0), color='purple', alpha=.3) + 
  geom_vline(xintercept=estimate) 

• The binomial distribution sampled by rbinom(m,k,1/2) is the distribution of the sum of \(k\) flips.
• Taking \(k=1\), we get the sum of one flip. A flip.
• Q. We’re getting more than 4 values. Why?

• This one combines the benefits of the last two.
  • We don’t call anyone twice.
  • The math is easy because flips are independent.
• But there’s a price. What is it?

• When we do this, our sample size is usually roughly half our population size.
• What should we do if we want to get a sample roughly a sixth our population size? Or a third?
• Are you equipped to do that?

• What do you think about this? How do you analyze it?
• Do you have enough intuition to sketch the distribution? Even a caricature?

estimate = NaN                                       
sam.plot = pop.plot + geom_hline(yintercept=estimate)    
sam.plot          

• It’s hard to shuffle six cards by hand, so I’m going to give you 24.
• Ace-6 in all suits. Ace is 1.
• Shuffle all of them, then take out the diamonds without changing their order.
• That’s 1-6 shuffled. The other cards were just filler to help you shuffle.

# editor-code-line-numbers: 2
by.shuffling = function(pop, sampling.rate=1/2) {
  N = round(m * sampling.rate) 
  shuffled = sample(1:m, replace=FALSE) 
  sam = pop[shuffled[1:N],]
  estimator(sam)
}

estimator.samples =  pops |> map_vec(by.shuffling)

ggplot(data.frame(sample=estimator.samples, rep=(1:R)/R)) + 
  geom_bar(aes(x=sample, y=after_stat(prop)), alpha=.3) +
  geom_point(aes(x=sample, y=rep), position=position_jitter(w=.01,h=0), color='purple', alpha=.3) + 
  geom_vline(xintercept=estimate) 

• This is how you shuffle in R. Why does it work? What does that tell you?

estimate = NaN                                       
sam.plot = pop.plot + geom_hline(yintercept=estimate)    
sam.plot          

• You probably don’t have the right die for this. Here’s a trick.
  • Roll a die with enough sides. Too many is fine.
  • Keep rolling until you get a number ≤ your sample size.
• The result of your last roll is as if you’d rolled the right die.
• This is called rejection sampling.

• Q. Why are we talking about this? Have we talked about a two-step scheme?