Objective : Use R as a calculator & Physical Notes.
\[ \text{Lec. Formulation : } \]
\[ \hat{\text{V}}(\hat{\tau_y})=\tau_x^2 \hat{\text{V}}(r) \\ \hat{\text{V}}(r)= \text{FPC}*\frac{1}{\mu_x^2}\frac{s_r^2}{n}\\ s^2_r=\frac{1}{n-1}\sum_{i\in S}(y_i - rx_i)\\ \text{FPC}=\frac{N-n}{N-1}\approx\frac{N-n}{N}\approx 1-\frac{n}{N} \]
x <- c(.3, .5, .4, .9, .7, .2, .6, .5, .8, .4, .8, .6)
y <- c(6, 9, 7, 19, 15, 5,12, 9, 20, 9, 18, 13)
T_x <- 75
sm_y <- sum(y)
sm_x <- sum(x)
df <- data.frame(x=x, y=y)
n <- 12
N <- 250
r <- sm_y/sm_x
T_y <- r*T_x
round(T_y, 2) # estimation
## [1] 1589.55
FPC <- 1-n/N
sr2 <- (1/(n-1))*sum((y-r*x)^2)
mu_x <- T_x/N
V_r <- FPC*(1/mu_x^2)*(sr2/n)
V_Ty <- V_r * T_x^2
B <- 2*sqrt(V_Ty)
round(B, 2) # Bound
## [1] 186.32
ci_Ty <- T_y + c(-1, 1)*B
ci_Ty # conf intvl
## [1] 1403.235 1775.870
T_y2 <- mean(y)*N
round(T_y2, 2) # est
## [1] 2958.33
V_Ty2 <-(var(y)/n)*FPC*N^2
B2 <- 2*sqrt(V_Ty2)
round(B2, 2) # bound
## [1] 730.13
# plt
plot(x,y)
abline(a=0, b = r, col = "green")
ci_Ty2 <- T_y2 + c(-1, 1)*B2
ci_Ty2 # conf intvl
## [1] 2228.199 3688.468
tab <- data.frame(
Estimate = c("T_y", "V_Ty", "B","T_y2", "V_Ty2", "B2"),
Value = c(T_y, V_Ty, B, T_y2, V_Ty2, B2)
)
tab
Given the origin-intercept, with slope r model good fit, it is clear that our ratio est is very likely better. Further given the fact that our est. bounds are dramatically smaller with our ratio est.
Further notice :
\[ N\bar{y} =\tau_x*r \iff \bar{y}=r\mu_x \\ \text{ or,}\\ N\bar{y} =\tau_x*r \iff \bar{x}=\mu_x \]
proof :
\[ N\bar{y} =\tau_xr \\ N\bar{y} =(N\mu_x)r \\ \boxed{\bar{y}=\mu_xr} \\ \bar{y}=\mu_x \frac{\bar{y}}{\bar{x}} \\ 1=\mu_x \frac{1}{\bar{x}}\\ \boxed{\bar{x}=\mu_x} \]
rm(list=ls()) #del gloabal env obj
x <- c("25,100", "32,200", "29,600", "35,000", "34,400", "26,500", "28,700", "28,200", "34,600", "32,700", "31,500", "30,600", "27,700", "28,500")
x <- gsub(",", "", x) |> as.numeric() # rep "," with empty space, then make it a num-vec
y <- c(3800, 5100, 4200, 6200, 5800, 4100, 3900, 3600, 3800, 4100, 4500, 5100, 4200, 4000)
df <- data.frame(x=x, y=y)
df
plot(x,y)
Recall :
\[ \hat{\text{V}}(r)= \text{FPC}*\frac{1}{\mu_x^2}\frac{s_r^2}{n}\\ s^2_r=\frac{1}{n-1}\sum_{i\in S}(y_i - rx_i)\\ \text{FPC}=\frac{N-n}{N-1}\approx\frac{N-n}{N}\approx 1-\frac{n}{N} \]
n <- 14
N <- 150
r <- mean(y)/mean(x)
round(r, 3)
## [1] 0.147
FPC <- 1-n/N
Sr2 <- (1/(n-1))*sum((y-r*x)^2)
V_r <- (Sr2/n)*(1/mean(x)^2)*FPC #est. of mu_x
B <- 2*sqrt(V_r)
round(B, 4)
## [1] 0.0102
ci <- r + c(-1,1)*B
ci
## [1] 0.1364745 0.1569654
recall :
\[ \hat{\text{V}}(\hat{\tau_y})=\tau_x^2 \hat{\text{V}}(r) \\ \hat{\text{V}}(r)= \text{FPC}*\frac{1}{\mu_x^2}\frac{s_r^2}{n}\\ s^2_r=\frac{1}{n-1}\sum_{i\in S}(y_i - rx_i)\\ \text{FPC}=\frac{N-n}{N-1}\approx\frac{N-n}{N}\approx 1-\frac{n}{N} \]
rm(list=ls())
x <- c(550, 720, 1500, 1020, 620, 980, 928, 1200, 1350, 1750, 670, 729, 1530)
y <- c(610, 780, 1600, 1030, 600, 1050, 977, 1440, 1570, 2210, 980, 865, 1710)
df <- data.frame(x=x, y=y)
df
plot(x,y)
T_x <- 128200
n <- 13
N <- 123
FPC <- (N-n)/(N-1)
mu_x <- T_x/N
r <- sum(y)/sum(x)
T_y <- T_x * r
round(T_y, 2) # est
## [1] 145943.8
y_hat <- r*x
sr2 <- (1/(n-1))*(sum((y-y_hat)^2))
V_r <- (sr2/n)*FPC*(1/mu_x^2)
V_Ty <- V_r*T_x^2
B <- 2*sqrt(V_Ty)
round(B, 2)
## [1] 7383.74