Tengo una línea de mejor ajuste. Necesito puntos de datos que no cambien mi línea de mejor ajuste

Estoy dando una presentación sobre líneas de ajuste. Tengo una función lineal simple, $y=1x+b$ . Estoy tratando de obtener puntos de datos dispersos que puedo poner en un diagrama de dispersión que mantendrá mi línea de mejor ajuste en la misma ecuación.

Me encantaría aprender esta técnica en R o Excel, lo que sea más fácil.

Elija cualquiera $(x_i)$ siempre que al menos dos de ellos difieran. Establezca una intersección $\beta_0$ y una pendiente $\beta_1$ y defina

Este ajuste es perfecto. Sin cambiar el ajuste, puede modificar $y_0$ a $y = y_0 + \varepsilon$ agregando cualquier vector de error $\varepsilon=(\varepsilon_i)$ siempre que sea ortogonal tanto al vector $x = (x_i)$ como al vector constante $(1,1,\ldots, 1)$ . Una manera fácil de obtener este tipo de error es escoger cualquier vector de $e$ y dejar que $\varepsilon$ sea los residuos sobre la regresión de $e$contra $x$ . En el siguiente código, $e$ se genera como un conjunto de valores normales aleatorios independientes con media $0$ y desviación estándar común.

Además, incluso puede preseleccionar la cantidad de dispersión, tal vez estipulando lo que debe ser $R^2$ . Dejando $\tau^2 = \text{var}(y_i) = \beta_1^2 \text{var}(x_i)$ , reescalar esos residuos para tener una varianza de

Este método es completamente general: todos los ejemplos posibles (para un conjunto dado de $x_i$ ) se pueden crear de esta manera.

Ejemplos

Cuarteto de Anscombe

Podemos reproducir fácilmente el Cuarteto de Anscombe de cuatro conjuntos de datos bivariados cualitativamente distintos que tienen las mismas estadísticas descriptivas (a través del segundo orden).

El código es notablemente simple y flexible.

set.seed(17)
rho <- 0.816                                             # Common correlation coefficient
x.0 <- 4:14
peak <- 10
n <- length(x.0)

# -- Describe a collection of datasets.
x <- list(x.0, x.0, x.0, c(rep(8, n-1), 19))             # x-values
e <- list(rnorm(n), -(x.0-peak)^2, 1:n==peak, rnorm(n))  # residual patterns
f <- function(x) 3 + x/2                                 # Common regression line

par(mfrow=c(2,2))
xlim <- range(as.vector(x))
ylim <- f(xlim + c(-2,2))
s <- sapply(1:4, function(i) {
  # -- Create data.
  y <- f(x[[i]])                                         # Model values
  sigma <- sqrt(var(y) * (1 / rho^2 - 1))                # Conditional S.D.
  y <- y + sigma * scale(residuals(lm(e[[i]] ~ x[[i]]))) # Observed values

  # -- Plot them and their OLS fit.
  plot(x[[i]], y, xlim=xlim, ylim=ylim, pch=16, col="Orange", xlab="x")
  abline(lm(y ~ x[[i]]), col="Blue")

  # -- Return some regression statistics.
  c(mean(x[[i]]), var(x[[i]]), mean(y), var(y), cor(x[[i]], y), coef(lm(y ~ x[[i]])))
})
# -- Tabulate the regression statistics from all the datasets.
rownames(s) <- c("Mean x", "Var x", "Mean y", "Var y", "Cor(x,y)", "Intercept", "Slope")
t(s)

La salida proporciona las estadísticas descriptivas de segundo orden para los datos $(x,y)$ de cada conjunto de datos. Las cuatro líneas son idénticas. Puede crear fácilmente más ejemplos alterando x(las coordenadas x) y e(los patrones de error) desde el principio.

Simulaciones

R$y$$\beta=(\beta_0,\beta_1)$$R^2$$0 \le R^2 \le 1$$x$

simulate <- function(x, beta, r.2) {
  sigma <- sqrt(var(x) * beta[2]^2 * (1/r.2 - 1))
  e <- residuals(lm(rnorm(length(x)) ~ x))
  return (y.0 <- beta[1] + beta[2]*x + sigma * scale(e))
}

(No sería difícil portar esto a Excel, pero es un poco doloroso).

$(x,y)$$60$ $x$$\beta=(1,-1/2)$$1$$-1/2$$R^2 = 0.5$

n <- 60
beta <- c(1,-1/2)
r.2 <- 0.5   # Between 0 and 1

set.seed(17)
x <- rnorm(n)

par(mfrow=c(1,4))
invisible(replicate(4, {
  y <- simulate(x, beta, r.2)
  fit <- lm(y ~ x)
  plot(x, y)
  abline(fit, lwd=2, col="Red")
}))

Al ejecutar summary(fit), puede verificar que los coeficientes estimados son exactamente como se especifica y los múltiples$R^2$es el valor pretendido Otras estadísticas, como el valor p de regresión, se pueden ajustar modificando los valores de$x_i$.