Mam linię najlepszego dopasowania. Potrzebuję punktów danych, które nie zmienią mojej linii najlepszego dopasowania

15

Prowadzę prezentację na temat dopasowania linii. Mam prostą funkcję liniową, . Próbuję uzyskać rozproszone punkty danych, które mogę umieścić na wykresie rozrzutu, który utrzyma moją linię najlepszego dopasowania na tym samym równaniu.y=1x+b

Chciałbym nauczyć się tej techniki w R lub Excel - w zależności od tego, co jest łatwiejsze.

Ryan Chase
źródło
1
Przypadek regresji wielokrotnej z dowolnym zestawem współczynników (z których twój jest przypadkiem szczególnym) omówiono w punkcie (2) tej odpowiedzi . Postępowanie zgodnie z instrukcjami rozwiązuje prosty przypadek regresji. Podejście to działa w prawie każdym pakiecie, w którym można symulować losowe wartości pożądanego rozkładu i dopasować modele regresji.
Glen_b
autodeskresearch.com/publications/samestats przedstawia niezłe uogólnienie: symulowane wyżarzanie służy do tworzenia wykresów rozrzutu, które nie tylko mają pożądane wartości statystyk podsumowujących, ale mają również określony kształt (np. „datasaurus”). Jest to praca Justina Matejki i George'a Fitzmaurice zatytułowana Same statystyki, różne wykresy: Generowanie zestawów danych o zróżnicowanym wyglądzie i statystyki identyczne poprzez symulowane wyżarzanie .
whuber

Odpowiedzi:

28

Wybierz dowolny (xi) pod warunkiem, że co najmniej dwa z nich się różnią. Ustaw punkt przecięcia β0 i nachylenie β1 i zdefiniuj

y0ja=β0+β1xja.

To dopasowanie jest idealne. Bez zmiany dopasowania można zmodyfikować y0 do y=y0+ε przez dodanie do niego dowolnego wektora błędu ε=(εja) , pod warunkiem że jest on ortogonalny zarówno do wektora x=(xja) i wektora stałego (1,1,,1) . Łatwym sposobem na uzyskanie takiego błędu jest wybranie dowolnego wektora mi i niech ε będzie resztą po regresji miprzeciw x . W poniższym kodzie mi jest generowany jako zbiór niezależnych losowych wartości normalnych ze średnią 0 i wspólnym odchyleniem standardowym.

Co więcej, można nawet preselekcji ilości rozproszonego, być może poprzez określenie co R2) powinno być. Niech τ2)=var(yja)=β12)var(xja) , przeskaluj te reszty, aby mieć wariancję

σ2)=τ2)(1/R2)-1).

Ta metoda jest w pełni ogólna: w ten sposób można utworzyć wszystkie możliwe przykłady (dla danego zestawu xja ).


Przykłady

Kwartet Anscombe

Możemy z łatwością odtworzyć Kwartet Anscombe czterech jakościowo odrębnych dwuwymiarowych zbiorów danych o takich samych statystykach opisowych (w drugim rzędzie).

Postać

Kod jest niezwykle prosty i elastyczny.

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)

Dane wyjściowe dają statystyki opisowe drugiego rzędu dla danych (x,y) dla każdego zestawu danych. Wszystkie cztery linie są identyczne. Możesz łatwo utworzyć więcej przykładów, zmieniając x(współrzędne x) i e(wzorce błędów) na początku.

Symulacje

Ryβ=(β0,β1)R2)0R2)1x

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))
}

(Przeniesienie tego do Excela nie byłoby trudne - ale jest to trochę bolesne).

(x,y)60 xβ=(1,-1/2))1-1/2)R2)=0,5

Postać

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")
}))

summary(fit)R2)xja

Whuber
źródło
1
Bardzo fajnie, dziękuję! Niestety, twoje podejście nie wydaje się mieć natychmiastowego zastosowania do tego pytania: zestawy danych podobne do Anscombe z tym samym polem i polem wąsów (średnia / standardowe / średnie / MAD / min / maks.) , Prawda?
Stephan Kolassa
@Stephan Masz rację, że tak nie jest, ponieważ jest to wysoce nieliniowy problem. Można go rozwiązać w podobny sposób - przede wszystkim poprzez znalezienie wykonalnych rozwiązań ograniczonego problemu optymalizacji - ale wymaga innej procedury optymalizacji, a rozwiązania nie są gwarantowane.
whuber