Jak symulować dane spełniające określone ograniczenia, takie jak posiadanie określonej średniej i odchylenia standardowego?

56

To pytanie jest motywowane moim pytaniem dotyczącym metaanalizy . Ale wyobrażam sobie, że przydałoby się to również w nauczaniu kontekstów, w których chcesz utworzyć zestaw danych, który dokładnie odzwierciedla istniejący opublikowany zestaw danych.

Wiem, jak generować losowe dane z danej dystrybucji. Na przykład, jeśli przeczytam o wynikach badania, które:

średnio 102,
odchylenie standardowe 5,2, oraz
wielkość próbki 72.

Mógłbym wygenerować podobne dane używając rnormR. Na przykład

set.seed(1234)
x <- rnorm(n=72, mean=102, sd=5.2)

Oczywiście średnia i SD nie byłyby dokładnie równe odpowiednio 102 i 5,2:

round(c(n=length(x), mean=mean(x), sd=sd(x)), 2)
##     n   mean     sd 
## 72.00 100.58   5.25

Ogólnie interesuje mnie, jak symulować dane, które spełniają zestaw ograniczeń. W powyższym przypadku stałymi są wielkość próby, średnia i odchylenie standardowe. W innych przypadkach mogą istnieć dodatkowe ograniczenia. Na przykład,

minimalne i maksymalne dane lub zmienna bazowa mogą być znane.
wiadomo, że zmienna przyjmuje tylko wartości całkowite lub tylko wartości nieujemne.
dane mogą obejmować wiele zmiennych o znanych wzajemnych korelacjach.

pytania

Ogólnie, jak mogę symulować dane, które dokładnie spełniają zestaw ograniczeń?
Czy są na ten temat artykuły? Czy są jakieś programy w R, które to robią?
Na przykład, w jaki sposób i powinienem symulować zmienną, aby miała określoną średnią i sd?

r dataset simulation random-generation Jeromy Anglim
źródło

1

Dlaczego chcesz, aby były dokładnie takie same jak opublikowane wyniki? Czy nie są to szacunki średniej populacji i odchylenia standardowego, biorąc pod uwagę ich próbkę danych. Biorąc pod uwagę niepewność tych szacunków, kto ma powiedzieć, że próbka, którą pokazałeś powyżej, nie jest zgodna z ich obserwacjami?

Gavin Simpson

4

Ponieważ wydaje się, że to pytanie zbiera odpowiedzi, które nie są zgodne ze znakiem (IMHO), chciałbym zauważyć, że koncepcyjnie odpowiedź jest prosta: ograniczenia równości są traktowane jak rozkłady krańcowe, a ograniczenia nierówności są wielowymiarowymi analogami obcięcia. Obcinanie jest stosunkowo łatwe w obsłudze (często z odrzucaniem próbkowania); trudniejszy problem polega na znalezieniu sposobu na próbkowanie tych rozkładów krańcowych. Oznacza to albo marginesy próbkowania, biorąc pod uwagę rozkład i ograniczenie, albo całkowanie, aby znaleźć rozkład krańcowy i próbkowanie na jego podstawie.

whuber

4

BTW, ostatnie pytanie jest banalne dla rodzin dystrybucji w skali lokalizacji. Np. x<-rnorm(72);x<-5.2*(x-mean(x))/sd(x)+102Robi lewę.

whuber

1

@ whuber, jak wspomina kardynał w komentarzu do mojej odpowiedzi (która wspomina o tej „sztuczce”) i komentarzu do innej odpowiedzi - ta metoda, ogólnie rzecz biorąc, nie zachowa zmiennych w tej samej rodzinie dystrybucyjnej, ponieważ dzielisz przez odchylenie standardowe próbki.

Makro,

5

@Macro To dobra uwaga, ale być może najlepszą odpowiedzią jest „oczywiście nie będą mieli takiej samej dystrybucji”! Rozkład, który chcesz, jest rozkładem zależnym od ograniczeń. Zasadniczo nie będzie to ta sama rodzina, co dystrybucja nadrzędna. Np. Każdy element próbki o wielkości 4 ze średnią 0 i SD 1 narysowany z rozkładu normalnego będzie miał prawie jednakowe prawdopodobieństwo na [-1,5, 1,5], ponieważ warunki ustalają górne i dolne granice możliwych wartości.

whuber

26

Ogólnie rzecz biorąc, aby średnia i wariancja w próbce były dokładnie równe z góry określonej wartości, można odpowiednio przesunąć i skalować zmienną. W szczególności, jeśli jest próbką, a następnie nowymi zmiennymi $X_1, X_2, ..., X_n$

Z_{i} = \sqrt{c_{1}} (\frac{X_{i} - \bar{X}}{s_{X}}) + c_{2}

$Z_i = \sqrt{c_{1}} \left( \frac{X_i-\overline{X}}{s_{X}} \right) + c_{2}$

gdzie $\overline{X} = \frac{1}{n} \sum_{i=1}^{n} X_i$ $s^{2}_{X} = \frac{1}{n-1} \sum_{i=1}^{n} (X_i - \overline{X})^2$ $Z_{i}$ $c_2$ $c_1$

B_{i} = a + (b - a) (\frac{X_{i} - min ({X_{1}, . . ., X_{n}})}{max ({X_{1}, . . ., X_{n}}) - min ({X_{1}, . . ., X_{n}})})

$B_i = a + (b-a) \left( \frac{ X_i - \min (\{X_1, ..., X_n\}) }{\max (\{X_1, ..., X_n\}) - \min (\{X_1, ..., X_n\}) } \right)$

utworzy zestaw danych który jest ograniczony do przedziału . $B_1, ..., B_n$ $(a,b)$

Uwaga: Te typy przesunięcia / skalowania zasadniczo zmienią rodzinę dystrybucyjną danych, nawet jeśli oryginalne dane pochodzą z rodziny o skali lokalizacji.

W kontekście tego rozkładu normalnegomvrnorm funkcja w R pozwala symulować normalne (lub wielowymiarowych normalny) dane z góry określonej próbki oznaczać / kowariancji przez ustawienie empirical=TRUE. W szczególności funkcja ta symuluje dane z rozkładu warunkowego zmiennej normalnie rozłożonej, biorąc pod uwagę średnią próbkę i (ko) wariancję równą z góry określonej wartości . Zauważ, że wynikowe rozkłady krańcowe nie są normalne, jak zauważył @whuber w komentarzu do głównego pytania.

Oto prosty przykład z jedną zmienną, w którym średnia próbki (z próbki ) jest ograniczona do 0, a odchylenie standardowe próbki wynosi 1. Widzimy, że pierwszy element jest znacznie bardziej podobny do rozkładu równomiernego niż normalny dystrybucja: $n=4$

library(MASS)
 z = rep(0,10000)
for(i in 1:10000)
{
    x = mvrnorm(n = 4, rep(0,1), 1, tol = 1e-6, empirical = TRUE)
    z[i] = x[1]
}
hist(z, col="blue")

$\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \$ wprowadź opis zdjęcia tutaj

Makro
źródło

1

nie będą normalnie rozmieszczone, że mogą one być w przybliżeniu tak, że rozmiar próbki jest duża. Pierwszy komentarz do odpowiedzi @ Seana odnosi się do tego.

Z_{i}

$Z_i$

kardynał

1

Cóż, to całkiem naturalna rzecz, którą chcesz zrobić ... i często nie sprawia to zbytniego problemu.

kardynał

1

+1. Nawiasem mówiąc , w tym przykładzie mundur jest dokładną odpowiedzią. (Widoczny spadek na końcach fabuły jest artefaktem rysowania histogramów przez R).

whuber

1

@ Whuber, dziękuję za motywowanie tego przykładu. Biorąc pod uwagę fakt, że rozkłady krańcowe zmieniają się, gdy uwarunkujesz średnią / wariancję próbki, wydaje się, że najlepszą „odpowiedzią” w duchu pytania PO jest po prostu symulacja danych ze średnią / wariancją populacji równą tej zgłoszonej jako próba ilości (jak sugeruje sam PO), prawda? W ten sposób otrzymujesz ilości próbek „podobne” do pożądanych, a rozkład krańcowy jest taki, jaki chciałeś.

Makro

1

@ whuber, Jeśli twoja próbka jest normalna, to ma rozkład , tak? Omawiana „nowa” zmienna będzie po prostu liniową kombinacją .

T_{i} = (X_{i} - \bar{X}) / s

$T_i = (X_i-\overline{X})/s$

t

$t$

T_{i}

$T_i$

Makro,

22

Jeśli chodzi o twoją prośbę o dokumenty, istnieją:

Chatterjee, S. & Firat, A. (2007). Generowanie danych z identycznymi statystykami, ale z odmienną grafiką: kontynuacja zbioru danych Anscombe. The American Statistician, 61 , 3, s. 248–254.

Nie jest to dokładnie to, czego szukasz, ale może służyć jako młyn do młyna.

Jest inna strategia, o której nikt chyba nie wspomniał. Możliwe jest generowanie (pseudo) losowych danych ze zbioru o rozmiarze tak że cały zestaw spełnia ograniczenia, o ile pozostałe danych jest ustawione na odpowiednie wartości. Wymagane wartości powinny być możliwe do rozwiązania za pomocą układu równań , algebry i pewnego smaru łokciowego. $N-k$ $N$ $k$ $k$ $k$

Na przykład, aby wygenerować zestaw danych z rozkładu normalnego, który będzie miał podaną średnią próbkową, i wariancję, , musisz naprawić wartości dwóch punktów: i . Ponieważ średnia próbki to: musi być: Przykładowa wariancja to: ten sposób (po zamianie powyższego na , foliowanie / dystrybucja i zmiana kolejności ... ) otrzymujemy: $N$ $\bar x$ $s^2$ $y$ $z$

\bar{x} = \frac{\sum_{i = 1}^{N - 2} x_{i} + y + z}{N}

$\bar x = \frac{\sum_{i=1}^{N-2}x_i\; + \;y\!+\!z}{N}$

y

$y$

y = N \bar{x} - (\sum_{i = 1}^{N - 2} x_{i} + z)

$y = N\bar x\; - \;\left(\sum_{i=1}^{N-2}x_i\!+\!z\right)$

s^{2} = \frac{\sum_{i = 1}^{N - 2} (x_{i} - \bar{x})^{2} + (y - \bar{x})^{2} + (z - \bar{x})^{2}}{N - 1}

$s^2 = \frac{\sum_{i=1}^{N-2}(x_i - \bar x)^2\; + \;(y - \bar x)^2\!+\!(z - \bar x)^2}{N-1}$

y

$y$

2 (N \bar{x} - \sum_{i = 1}^{N - 2} x_{i}) z - 2 z^{2} = N {\bar{x}}^{2} (N - 1) + \sum_{i = 1}^{N - 2} x_{i}^{2} + {[\sum_{i = 1}^{N - 2} x_{i}]}^{2} - 2 N \bar{x} \sum_{i = 1}^{N - 2} x_{i} - (N - 1) s^{2}

$2(N\bar{x}\! - \!\sum_{i=1}^{N-2}x_i)z - 2z^2 = N\bar{x}^2(N\!-\!1) + \sum_{i=1}^{N-2}x_i^2 + \left[\sum_{i=1}^{N-2}x_i\right]^2 - 2N\bar{x}\sum_{i=1}^{N-2}x_i - (N\!-\!1)s^2$ Jeśli weźmiemy , , oraz jako negacja RHS, możemy rozwiązać dla za pomocą wzoru kwadratowego . Na przykład można użyć następującego kodu:

a = - 2

$a=-2$

b = 2 (N \bar{x} - \sum_{i = 1}^{N - 2} x_{i})

$b=2(N\bar{x} - \sum_{i=1}^{N-2}x_i)$

c

$c$

z

$z$ R

find.yz = function(x, xbar, s2){
  N    = length(x) + 2
  sumx = sum(x)
  sx2  = as.numeric(x%*%x)          # this is the sum of x^2
  a    = -2
  b    = 2*(N*xbar - sumx)
  c    = -N*xbar^2*(N-1) - sx2 - sumx^2 + 2*N*xbar*sumx + (N-1)*s2
  rt   = sqrt(b^2 - 4*a*c)

  z    = (-b + rt)/(2*a)
  y    = N*xbar - (sumx + z)
  newx = c(x, y, z)
  return(newx)
}

set.seed(62)
x    = rnorm(2)
newx = find.yz(x, xbar=0, s2=1)
newx                                # [1] 0.8012701  0.2844567  0.3757358 -1.4614627
mean(newx)                          # [1] 0
var(newx)                           # [1] 1

Podejście to należy zrozumieć. Po pierwsze, nie ma gwarancji, że zadziała. Na przykład, możliwe jest, że początkowe Dane są takie, że nie ma wartości i tego, że istnieją będzie wariancji Otrzymaną równy . Rozważać: $N-2$ $y$ $z$ $s^2$

set.seed(22)    
x    = rnorm(2)
newx = find.yz(x, xbar=0, s2=1)
Warning message:
In sqrt(b^2 - 4 * a * c) : NaNs produced
newx                                # [1] -0.5121391  2.4851837        NaN        NaN
var(c(x, mean(x), mean(x)))         # [1] 1.497324

Po drugie, podczas gdy standaryzacja sprawia, że rozkłady krańcowe wszystkich twoich zmiennych są bardziej jednolite, to podejście wpływa tylko na dwie ostatnie wartości, ale sprawia, że ich rozkłady krańcowe są wypaczone:

set.seed(82)
xScaled = matrix(NA, ncol=4, nrow=10000)
for(i in 1:10000){
  x           = rnorm(4)
  xScaled[i,] = scale(x)
}

(wstaw działkę)

set.seed(82)
xDf = matrix(NA, ncol=4, nrow=10000)
i   = 1
while(i<10001){
  x       = rnorm(2)
  xDf[i,] = try(find.yz(x, xbar=0, s2=2), silent=TRUE)  # keeps the code from crashing
  if(!is.nan(xDf[i,4])){ i = i+1 }                      # increments if worked
}

(wstaw działkę)

Po trzecie, uzyskana próbka może nie wyglądać bardzo normalnie; może to wyglądać tak, jakby zawierało „wartości odstające” (tj. punkty, które pochodzą z innego procesu generowania danych niż reszta), ponieważ tak jest w istocie. Jest to mniej prawdopodobne, że będzie to stanowić problem w przypadku większych próbek, ponieważ statystyki próbek z wygenerowanych danych powinny być zbieżne z wymaganymi wartościami, a zatem wymagać mniejszego dostosowania. Przy mniejszych próbkach zawsze możesz połączyć to podejście z algorytmem akceptowania / odrzucania, który próbuje ponownie, jeśli wygenerowana próbka ma statystyki kształtu (np. Skośność i kurtoza), które są poza dopuszczalnymi granicami (por. Komentarz kardynała ) lub rozszerzyć takie podejście do generowania próbki o ustalonej średniej, wariancji, skośności ikurtoza (ale zostawię tobie algebrę). Alternatywnie, możesz wygenerować niewielką liczbę próbek i użyć tej z najmniejszą (powiedzmy) statystyką Kołmogorowa-Smirnowa.

library(moments)
set.seed(7900)  
x = rnorm(18)
newx.ss7900 = find.yz(x, xbar=0, s2=1)
skewness(newx.ss7900)                       # [1] 1.832733
kurtosis(newx.ss7900) - 3                   # [1] 4.334414
ks.test(newx.ss7900, "pnorm")$statistic     # 0.1934226

set.seed(200)  
x = rnorm(18)
newx.ss200 = find.yz(x, xbar=0, s2=1)
skewness(newx.ss200)                        # [1] 0.137446
kurtosis(newx.ss200) - 3                    # [1] 0.1148834
ks.test(newx.ss200, "pnorm")$statistic      # 0.1326304 

set.seed(4700)  
x = rnorm(18)
newx.ss4700 = find.yz(x, xbar=0, s2=1)
skewness(newx.ss4700)                       # [1]  0.3258491
kurtosis(newx.ss4700) - 3                   # [1] -0.02997377
ks.test(newx.ss4700, "pnorm")$statistic     # 0.07707929S

(dodaj działkę)

gung - Przywróć Monikę
źródło

10

Ogólna technika to „Metoda odrzucenia”, w której po prostu odrzucasz wyniki, które nie spełniają twoich ograniczeń. Jeśli nie masz wskazówek (np. MCMC), możesz generować wiele przypadków (w zależności od scenariusza), które są odrzucane!

Jeśli szukasz czegoś w rodzaju średniej i odchylenia standardowego i możesz utworzyć metrykę odległości, aby określić, jak daleko jesteś od celu, możesz użyć optymalizacji do wyszukiwania zmiennych wejściowych, które dają pożądany wynik wartości.

Jako brzydki przykład, w którym szukamy losowego jednorodnego wektora o długości 100, który ma średnią = 0 i odchylenie standardowe = 1.

# simplistic optimisation example
# I am looking for a mean of zero and a standard deviation of one
# but starting from a plain uniform(0,1) distribution :-)
# create a function to optimise
fun <- function(xvec, N=100) {
  xmin <- xvec[1]
  xmax <- xvec[2]
  x <- runif(N, xmin, xmax)
  xdist <- (mean(x) - 0)^2 + (sd(x) - 1)^2
  xdist
}
xr <- optim(c(0,1), fun)

# now lets test those results
X <- runif(100, xr$par[1], xr$par[2])
mean(X) # approx 0
sd(X)   # approx 1

Sean
źródło

7

Ograniczenia występujące z prawdopodobieństwem zerowym są trudne do spełnienia. ;-) W przypadku konkretnego przykładu odpowiednia zmiana i dylatacja z łatwością osiągają wyznaczone cele, choć można chcieć przeanalizować nieco głębiej, aby zobaczyć, jak taka operacja zakłóca rozkład danych.

kardynał

Dzięki. Z pewnością łatwo byłoby odrzucić obserwacje mniejsze niż min i większe niż maks. Widzę, jak można zdefiniować to jako problem optymalizacji. Byłoby wspaniale zobaczyć kilka przykładów, a może sugestie co dalej.

Jeromy Anglim

1

@cardinal - uzgodniony. Należy spojrzeć na rozkłady (tj. Histogram) zarówno wejściowych liczb symulowanych, jak i wyjściowych, ponieważ czasami mogą one wyglądać naprawdę bardzo dziwnie!

Sean

9

Czy są jakieś programy w R, które to robią?

Pakiet Runuran R zawiera wiele metod generowania losowych zmiennych. Wykorzystuje biblioteki C z projektu UNU.RAN (Universal Non-Uniform RAndom Number generator). Moja własna wiedza na temat generowania zmiennych losowych jest ograniczona, ale winieta Runuran zapewnia ładny przegląd. Poniżej znajdują się dostępne metody w pakiecie Runuran, zaczerpnięte z winiety:

Ciągłe dystrybucje:

Adaptacyjne odrzucanie próbek
Odwrotne odwrócenie transformowanej gęstości
Wielomianowa interpolacja odwrotnego CDF
Prosta metoda stosunku mundurów
Odrzucone przekształcenie gęstości

Dyskretne rozkłady:

Dyskretna automatyczna inwersja odrzucania
Metoda Alias-Urn
Metoda tabeli dyskretnej inwersji

Rozkłady wielowymiarowe:

Algorytm Hit-and-Run z metodą Ratio-of-Uniforms
Metoda wielowymiarowego naiwnego stosunku mundurów

Przykład:

Na przykład, załóżmy, że chcesz wygenerować rozkład normalny ograniczony od 0 do 100:

require("Runuran")

## Normal distribution bounded between 0 and 100
d1 <- urnorm(n = 1000, mean = 50, sd = 25, lb = 0, ub = 100)

summary(d1)
sd(d1)
hist(d1)

Ta urnorm()funkcja jest wygodną funkcją owijania. Uważam, że za kulisami używa metody interpolacji wielomianowej odwrotnej CDF, ale nie jestem pewien. Dla czegoś bardziej złożonego, powiedzmy, dyskretny rozkład normalny ograniczony od 0 do 100:

require("Runuran")

## Discrete normal distribution bounded between 0 and 100
# Create UNU.RAN discrete distribution object
discrete <- unuran.discr.new(pv = dnorm(0:100, mean = 50, sd = 25), lb = 0, ub = 100)

# Create UNU.RAN object using the Guide-Table Method for Discrete Inversion
unr <- unuran.new(distr = discrete, method = "dgt")

# Generate random variates from the UNU.RAN object
d2 <- ur(unr = unr, n = 1000)

summary(d2)
sd(d2)
head(d2)
hist(d2)

jthetzel
źródło

3

Wygląda na to, że pakiet R spełniający twoje wymagania został opublikowany wczoraj! simstudy Keith Goldfeld

Symuluje zestawy danych w celu poznania technik modelowania lub lepszego zrozumienia procesów generowania danych. Użytkownik określa zestaw relacji między zmiennymi towarzyszącymi i generuje dane na podstawie tych specyfikacji. Ostateczne zestawy danych mogą reprezentować dane z randomizowanych prób kontrolnych, powtarzalnych pomiarów (podłużnych) i losowych prób skupień. Brak można wygenerować za pomocą różnych mechanizmów (MCAR, MAR, NMAR).

Tyelcie
źródło

1

Ani w winiecie, ani na stronie głównej programu nie wymieniono dokładnego spełnienia ograniczeń. Jak myślisz, dlaczego ten pakiet spełnia wymóg czerpania z dystrybucji warunkowych?

gg

2

To odpowiedź przychodzi tak późno, że prawdopodobnie nie ma ona znaczenia, ale zawsze istnieje rozwiązanie MCMC. Mianowicie, aby rzutować gęstość połączenia próbki na kolektor zdefiniowany przez ograniczenia, na przykład Jedynym problemem jest wówczas symulacja wartości w tym kolektorze, tj. znalezienie parametryzacji właściwego wymiaru. Artykuł z 2015 roku autorstwa Bornna, Shepharda i Solgi analizuje ten bardzo problem (z ciekawą, jeśli nie ostateczną odpowiedzią ).

\prod_{i = 1}^{n} f (x_{i})

$\prod_{i=1}^n f(x_i)$

\sum_{i = 1}^{n} x_{i} = μ_{0} \sum_{i = 1}^{n} x_{i}^{2} = σ_{0}^{2}

$\sum_{i=1}^n x_i=\mu_0\qquad\sum_{i=1}^n x_i^2=\sigma_0^2$

Xi'an
źródło

2

Ta odpowiedź rozważa inne podejście do przypadku, w którym chcesz zmusić wariaty do położenia się w określonym zakresie i dodatkowo dyktować średnią i / lub wariancję.

Ogranicz naszą uwagę do interwału jednostkowego . Użyjmy średniej ważonej dla ogólności, więc napraw niektóre wagi pomocą lub ustaw jeśli chcesz standardowej wagi. Załóżmy, że ilości i reprezentują odpowiednio pożądaną (ważoną) średnią i (ważoną) wariancję. Górna granica jest konieczna, ponieważ jest to maksymalna możliwa wariancja w jednostkowym przedziale. Jesteśmy zainteresowani rysowaniem niektórych wariantów z z tymi ograniczeniami momentu. $[0,1]$ $w_k\in[0,1]$ $\sum_{k=1}^Nw_k=1$ $w_k=1/N$ $\mu\in(0,1)$ $0<\sigma^2<\mu(1-\mu)$ $\sigma^2$ $x_1,...,x_N$ $[0,1]$

Najpierw narysujemy kilka odmian z dowolnego rozkładu, np. . Rozkład ten wpłynie na kształt rozkładu końcowego. Następnie ograniczamy je do przedziału jednostek za pomocą funkcji logistycznej: $y_1,...,y_N$ $N(0,1)$ $[0,1]$

x_{k} = \frac{1}{1 + e^{- (y_{k} v - h)}}

$x_k=\frac{1}{1+e^{-(y_k v-h)}}$

Jednak zanim to zrobimy, jak widać w powyższym równaniu, przekształcamy za pomocą translacji i skali . Jest to analogiczne do pierwszego równania w odpowiedzi @ Macro. Sztuką jest wybrać i tak, że przekształcone zmienne mieć pożądaną chwilę (y). Oznacza to, że do przechowywania potrzebujemy jednego lub obu następujących elementów: $y_k$ $h$ $v$ $h$ $v$ $x_1,...,x_N$

μ = \sum_{k = 1}^{N} \frac{w_{k}}{1 + e^{- (y_{k} v - h)}} σ^{2} = \sum_{k = 1}^{N} \frac{w_{k}}{(1 + e^{- (y_{k} v - h)})^{2}} - {(\sum_{k = 1}^{N} \frac{w_{k}}{1 + e^{- (y_{k} v - h)}})}^{2}

$\mu=\sum_{k=1}^N \frac{w_k}{1+e^{-(y_k v-h)}} \\ \sigma^2=\sum_{k=1}^N \frac{w_k}{(1+e^{-(y_k v-h)})^2} - \left( \sum_{k=1}^N \frac{w_k}{1+e^{-(y_k v-h)}} \right)^2$

Analityczne odwrócenie tych równań dla i nie jest możliwe, ale wykonanie liczbowe jest proste, zwłaszcza że pochodne w odniesieniu do i są łatwe do obliczenia; zajmuje tylko kilka iteracji metody Newtona. $v$ $h$ $v$ $h$

Jako pierwszy przykład załóżmy, że zależy nam jedynie na ograniczeniu średniej ważonej, a nie na wariancji. Fix , , , . Następnie dla leżących u podstaw rozkładów , i otrzymujemy odpowiednio następujące histogramy i takie, że średnia zmiennych wynosi dokładnie (nawet dla małych ): $\mu=0.8$ $v=1$ $w_k=1/N$ $N=200000$ $N(0,1)$ $N(0,0.1)$ $\text{Unif}(0,1)$ $0.8$ $N$

Następnie ograniczmy zarówno średnią, jak i wariancję. Weźmy , , i rozważmy trzy pożądane odchylenia standardowe . Stosując ten sam podstawowy rozkład , tutaj są histogramy dla każdego: $\mu=0.2$ $w_k=1/N$ $N=2000$ $\sigma=0.1,0.05,0.01$ $N(0,1)$

Pamiętaj, że mogą one wyglądać nieco rozproszone w wersji beta, ale nie są.

Ian Hincks
źródło

1

W mojej odpowiedzi tutaj wymieniłem trzy pakiety R do zrobienia tego:

abalter
źródło

Link do referencji musi zawierać format. Czy powinien to być komentarz?

abalter

Jak symulować dane spełniające określone ograniczenia, takie jak posiadanie określonej średniej i odchylenia standardowego?

pytania

Odpowiedzi: