Wenn ich Ihre Frage richtig verstehe, ist das ganz einfach. Sie müssen nur entscheiden, welche Verteilung Ihre Fehler haben sollen, und die entsprechende Zufallsgenerierungsfunktion verwenden.
Es gibt eine Reihe von verzerrten Verteilungen, daher müssen Sie herausfinden, welche Sie mögen. Darüber hinaus sind die meisten verzerrten Verteilungen (z. B. logarithmisch normal, Chi-Quadrat, Gamma, Weibull usw.) recht verzerrt, sodass einige geringfügige Anpassungen erforderlich wären (z. B. mit multiplizieren ). - 1
Hier ist ein Beispiel zum Ändern Ihres Codes:
set.seed(5840) # this makes the example exactly reproducible
N <- 100
x <- rnorm(N)
beta <- 0.4
errors <- rlnorm(N, meanlog=0, sdlog=1)
errors <- -1*errors # this makes them left skewed
errors <- errors - 1 # this centers the error distribution on 0
y <- 1 + x*beta + errors
Ich sollte an dieser Stelle beachten, dass die Regression keine Annahmen über die Verteilungen von oder , sondern nur über die Fehler (siehe hier: Was ist, wenn die Residuen normal verteilt sind, y jedoch nicht? ). Das war also der Schwerpunkt meiner obigen Antwort. Y εX.Y.ε
Update: Hier ist eine rechtwinklige Version mit den als Weibull verteilten Fehlern:
set.seed(5840) # this makes the example exactly reproducible
N <- 100
x <- rnorm(N)
beta <- 0.4
errors <- rweibull(N, shape=1.5, scale=1)
# errors <- -1*errors # this makes them left skewed
errors <- errors - factorial(1/1.5) # this centers the error distribution on 0
y <- 1 + x*beta + errors
Weibull-Daten sind bereits richtig verzerrt, sodass wir ihre Richtung nicht ändern müssen (dh wir lassen das -1*errors
Teil fallen). Auf der Wikipedia-Seite für die Weibull-Verteilung sehen wir auch, dass der Mittelwert eines Weibull sein sollte:. Wir wollen diesen Wert von jedem der Fehler subtrahieren, damit die resultierende Fehlerverteilung auf zentriert ist . Dadurch kann der strukturelle Teil (dh ) Ihres Codes den strukturellen Teil des Datengenerierungsprozesses genau wiedergeben. 0E.[ W.] = ( 1 / s h a p e ) !01 + x*beta
Die ExGaußsche Verteilung ist die Summe aus Normal und Exponential. Es gibt eine Funktion ? RexGAUS im Paket gamlss.dist , um diese zu generieren. Ich habe dieses Paket nicht, aber Sie sollten in der Lage sein, meinen obigen Code ohne allzu große Schwierigkeiten anzupassen. Sie können auch eine zufällige Normalvariable (via rnorm()
) und eine Exponentialvariable (via rexp()
) generieren und diese ganz einfach summieren. Denken Sie daran, den Populationsmittelwert von jedem Fehler zu subtrahieren, bevor Sie die Fehler zum strukturellen Teil des Datenerzeugungsprozesses hinzufügen. (Achten Sie darauf , nicht die subtrahieren Probe Mittelwert , aber!) μ+1/λmean(errors)
Einige abschließende, nicht verwandte Kommentare: Ihr Beispielcode in der Frage ist etwas durcheinander (was keine Beleidigung bedeutet). Da rnorm(N)
Daten mit mean=0
und sd=1
standardmäßig 0.4*rnorm(N)
generiert , werden generiert rnorm(N, mean=0, sd=0.4)
. Ihr Code (und möglicherweise Ihr Denken) wird viel klarer, wenn Sie die letztere Formulierung verwenden. Außerdem beta
scheint Ihr Code für verwirrt zu sein. Wir denken im Allgemeinen an dieβin einem Modell vom Regressionstyp als Parameter, nicht als Zufallsvariable. Das heißt, es ist eine unbekannte Konstante, die das Verhalten des Datenerzeugungsprozesses bestimmt, aber die stochastische Natur des Prozesses wird durch die Fehler eingekapselt. Dies ist nicht die Art und Weise, wie wir darüber denken, wenn wir mit Mehrebenenmodellen arbeiten, und Ihr Code scheint auf halbem Weg zwischen einem Standardregressionsmodell und dem Code für ein Mehrebenenregressionsmodell zu liegen. Die separate Angabe Ihrer Betas ist eine gute Idee, um die konzeptionelle Klarheit des Codes zu erhalten. Bei einem Standard-Regressionsmodell würden Sie jedoch jeder Beta (z beta0 <- 1; beta1 <- .04
. B. ) nur eine einzige Nummer zuweisen .