Symulowanie wycieczki Browna przy użyciu mostu Browna?

11

Chciałbym zasymulować proces wycieczki Browna (ruch Browna, który jest warunkowany, zawsze jest dodatni, gdy do przy ). Ponieważ proces wycieczki Browna jest mostem Browna, który jest uwarunkowany, aby zawsze był dodatni, miałem nadzieję symulować ruch wycieczki Browna za pomocą mostu Browna.0 t = 10<t<10t=1

W R używam pakietu „e1017” do symulacji procesu mostu Browna. Jak mogę użyć tego procesu mostu Browna, aby utworzyć wycieczkę Browna?

RPz
źródło
4
Czy nie wystarczy zasymulować bezwzględną wartość mostu Browna?
Alex R.
1
@AlexR. no [padding]
P.Windridge
1
Warto jednak przypomnieć, że ruch Browna uwarunkowany dodatnim może być zrealizowany poprzez odzwierciedlenie BM wokół jego maksymalnego działania, co jest wynikiem Pitmana. Innym sposobem na osiągnięcie BM uwarunkowanego pozostaniem dodatnim jest bezwzględna wartość 3d BM .
P.Windridge
1
@AlexR. - Zaktualizowałem moją odpowiedź poniżej, aby pokazać, że nawet w przypadku prostych, losowych spacerów, pozytywne wyniki warunkowania wywołują różne zachowania niż przyjmowanie wartości bezwzględnej. W przypadku mostów Browna, intuicyjnie dla małych , zachowanie B B t jest jak | W t (ponieważ B B t = W t - t W 1 ) i BM spełniają prawo iterowanego logarytmu (więc „ O p ( t ) ” jest nieistotne dla wystarczająco małego t . Zatem | B B ttBBt|WtBBt=WttW1Op(t)tjest jak odbity BM dla małego t . To ma zupełnie inne zachowanie, aby W t warunkowane pozostają pozytywne ...|BBt|tWt
P.Windridge

Odpowiedzi:

7

Wycieczka Browna może być zbudowana z mostu przy użyciu następującej konstrukcji Vervaat: https://projecteuclid.org/download/pdf_1/euclid.aop/1176995155

Szybkie przybliżenie w R przy użyciu kodu BB @ whubera to

n <- 1001
times <- seq(0, 1, length.out=n)

set.seed(17)
dW <- rnorm(n)/sqrt(n)
W <- cumsum(dW)

# plot(times,W,type="l") # original BM

B <- W - times * W[n]   # The Brownian bridge from (0,0) to (1,target)

# plot(times,B,type="l")

# Vervaat construction
Bmin <- min(B)
tmin <- which(B == Bmin)
newtimes <- (times[tmin] + times) %% 1
J<-floor(newtimes * n)
BE <- B[J] - Bmin
plot(1:length(BE)/n,BE,type="l")

wprowadź opis zdjęcia tutaj

Oto kolejna fabuła (z set.seed (21)). Kluczową obserwacją podczas wycieczki jest to, że uwarunkowanie faktycznie manifestuje się jako „odpychanie” od 0, i jest mało prawdopodobne, aby wycieczka zbliżyła się do wewnątrz ( 0 , 1 ) . 0(0,1)wprowadź opis zdjęcia tutaj


Poza tym: rozkład wartości bezwzględnej mostu Browna i wycieczka, ( B B t ) 0 t 1 uwarunkowane jako dodatnie , nie są takie same. Intuicyjnie wycieczka jest odpychana od źródła, ponieważ ścieżki Browna, które znajdują się zbyt blisko źródła, prawdopodobnie wkrótce staną się ujemne, a zatem będą karane przez warunkowanie.(|BBt|)0t1(BBt)0t1

Można to zilustrować prostym pomostem i wycieczką na stopniach, który jest naturalnym dyskretnym analogiem BM (i zbiega się w BM, gdy kroki stają się duże i przeskalujesz).6

Rzeczywiście, weź symetryczny SRW, zaczynając od . Najpierw zastanówmy się nad warunkowaniem „pomostowym” i zobaczmy, co się stanie, jeśli weźmiemy wartość bezwzględną. Rozważyć wszystkie prostych ścieżek s o długości 6 , które zaczynają się i kończą w 0 . Liczba takich ścieżek wynosi . Są z nich, dla których . Innymi słowy, prawdopodobieństwo, że wartość bezwzględna naszego „mostka” SRW (uwarunkowanego zakończeniem na ) będzie miała wartość 0 na etapie wynosi .0s602× ( 4(63)=20| s2| =00212/20=0,62×(42)=12|s2|=00212/20=0.6

Po drugie, rozważymy uwarunkowanie „wycieczki”. Liczba nieujemna prostych odcinków o długości , które kończą się w jest liczbą kataloński . Dokładnie z tych ścieżek mają . Zatem prawdopodobieństwo, że nasza „wycieczka” SRW (uwarunkowana pozostanie dodatnia i zakończy się na ) o wartości 0 w kroku wynosi .6 = 2 3 0 C m = 3 = ( 2 ms6=2302s2=0022/5=0,4<0,6Cm=3=(2mm)/(m+1)=52s2=0022/5=0.4<0.6

Jeśli nadal masz wątpliwości, że zjawisko to utrzymuje się w granicach, możesz rozważyć prawdopodobieństwo dla mostów SRW i skoków długości uderzających 0 w kroku .2 n4n2n

Na wycieczkę SRW: mamy przy użyciu aysmptotics z wikipedii https://en.wikipedia.org/wiki / Catalan_number . To znaczy, że w końcu jest jak .Cn - 3 / 2

P(S2n=0|Sj0,j4n,S4n=0)=Cn2/C2n(42n/πn3)/(42n/(2n)3π)
cn3/2

Dla abs (most SRW): przy użyciu asymptotyków z wikipedii https://en.wikipedia.org/wiki/Binomial_coefficient . To jest jak .Cn - 1 / 2

P(|S2n|=0|S4n=0)=(2nn)2/(4n2n)(4n/πn)2/(42n/2nπ)
cn1/2

Innymi słowy, asymptotyczne prawdopodobieństwo zobaczenia warunku SRW dodatniego przy pobliżu środka jest znacznie mniejsze niż dla wartości bezwzględnej mostu. 0


Oto alternatywna konstrukcja oparta na procesie 3D Bessela zamiast mostu Browna. Korzystam z faktów wyjaśnionych w https://projecteuclid.org/download/pdf_1/euclid.ejp/1457125524

Omówienie - 1) Symuluj proces 3d Bessela. To jest jak BM uwarunkowane byciem pozytywnym. 2) Zastosuj odpowiednie przeskalowanie czasoprzestrzeni w celu uzyskania mostu Bessela 3 (równanie (2) w pracy). 3) Wykorzystaj fakt (zauważony tuż po Twierdzeniu 1 w pracy), że most Bessel 3 ma faktycznie taki sam rozkład jak wycieczka Browna.

Nieznaczną wadą jest to, że trzeba uruchomić proces Bessela przez dość długi czas (T = 100 poniżej) na stosunkowo drobnej siatce, aby skalowanie czasoprzestrzeni zaczęło się na końcu.

## Another construction of Brownian excursion via Bessel processes
set.seed(27092017)
## The Bessel process must run for a long time in order to construct a bridge
T <- 100
n <- 100001
d<-3 # dimension for Bessel process
dW <- matrix(ncol = n, nrow = d, data=rnorm(d*n)/sqrt(n/T))
dW[,1] <- 0
W <- apply(dW, 1, cumsum)
BessD <- apply(W,1,function(x) {sqrt(sum(x^2))})

times <- seq(0, T, length.out=n)
# plot(times,BessD, type="l") # Bessel D process


times01 <- times[times < 1]
rescaletimes <- pmin(times01/(1-times01),T)
# plot(times01,rescaletimes,type="l") # compare rescaled times

# create new time index
rescaletimeindex <- sapply(rescaletimes,function(x){max(which(times<=x))} )

BE <- (1 - times01) * BessD[rescaletimeindex]
plot(times01,BE, type="l")

Oto wynik: wprowadź opis zdjęcia tutaj

P.Windridge
źródło
5

Odbicie Zasada twierdzi

jeżeli ścieżka procesu Wienera osiąga wartość w czasie , to kolejna ścieżka po czasie ma taki sam rozkład, jak odbicie kolejnej ścieżki o wartościf(t)f(s)=at=ssa

Wikipedia , dostęp 26.09.2017.

W związku z tym możemy symulować most Browna i odzwierciedlić go o wartości po prostu przyjmując jego wartość bezwzględną. Most Browna jest symulowany przez odjęcie trendu od punktu początkowego do końca od samego ruchu Browna . (Bez utraty ogólności możemy mierzyć czas w jednostkach, które sprawiają, że Zatem w czasie po prostu odejmij od .)a=0(0,0)(T,B(T))BT=1tB(T)tB(t)

Tę samą procedurę można zastosować, aby wyświetlić ruch Browna pod warunkiem nie tylko powrotu do określonej wartości w czasie (wartość wynosi dla mostka), ale także pozostania między dwiema granicami (które koniecznie obejmują wartość początkową od w czasie i określonej wartości końcowej).T>0000

![Postać

Ten ruch Browna zaczyna się i kończy na wartości zero: jest to Most Browna.

Rysunek 2

Czerwony wykres to wycieczka Browna opracowana na podstawie poprzedniego mostu Browna: wszystkie jej wartości są nieujemne. Niebieski wykres opracowano w ten sam sposób, odzwierciedlając most Browna między kropkowanymi liniami za każdym razem, gdy je napotyka. Szary wykres pokazuje oryginalny most Browna.

Obliczenia są proste i szybkie: podziel zestaw czasów na małe interwały, wygeneruj niezależne, identycznie rozmieszczone normalne przyrosty dla każdego interwału, kumuluj je, odejmuj trend i wykonuj potrzebne refleksje.

Oto Rkod. W nim Wjest oryginalny ruch Browna, Bjest most Browna i B2jest ograniczony między dwiema określonymi wartościami ymin(nie dodatnimi) i ymax(nieujemnymi). Jego technika wykonywania odbicia przy użyciu %%operatora modułu i minimum komponentowego pminmoże być praktyczna.

#
# Brownian bridge in n steps from t=0 to t=1.
#
n <- 1001
times <- seq(0, 1, length.out=n)
target <- 0                        # Constraint at time=1
set.seed(17)
dW <- rnorm(n)
W <- cumsum(dW)
B <- W + times * (target - W[n])   # The Brownian bridge from (0,0) to (1,target)
#
# The constrained excursion.
#
ymax <- max(abs(B))/5              # A nice limit for illustration
ymin <- -ymax * 2                  # Another nice limit
yrange2 <- 2*(ymax - ymin)
B2 <- (B - ymin) %% yrange2
B2 <- pmin(B2, yrange2-B2) + ymin
Whuber
źródło
czy mógłbyś udostępnić kod swojej „wycieczki po Brownie” (czerwony wykres) Na oko wygląda bardziej jak rodzaj odbijanego ruchu Browna ograniczonego do . Myślę, że ma to nieco inny rozkład niż wycieczka, która doświadcza odpychania od pochodzenia, tzn. Wasza realizacja (na czerwono) wydaje się raczej nietypowa. 0
P.Windridge
@ P.Windridge Przepraszamy; Zapomniałem: wycieczka jest abs(B). Pamiętaj, że ma to być ruch Browna uwarunkowany dwoma ograniczeniami: jest równy targetczasowi i wszędzie jest nieujemny. 1
whuber
1
Nie zapomniałem :) Mówię, że uważam, że ma raczej inny rozkład niż uwarunkowany do bądź pozytywny (tj. wycieczka) :) ( B B t ) 0 t 1(abs(BBt))0t1(BBt)0t1
P.Windridge
4
Dystrybucje są różne, więc głosuję w dół przepraszam.
P.Windridge
2

Możesz użyć metody odrzucenia: symuluj mosty Browna i zachowaj te pozytywne. To działa.

Ale. Jest bardzo powolny, ponieważ wiele przykładowych trajektorii jest odrzucanych. Im większa ustawiona „częstotliwość”, tym mniejsze prawdopodobieństwo znalezienia trajektorii.

succeeded <- FALSE
while(!succeeded)
{
  bridge <- rbridge(end = 1, frequency = 500)
  succeeded=all(bridge>=0)
}
plot(bridge)

Możesz to przyspieszyć, utrzymując również negatywne trajektorie.

while(!succeeded)
{
  bridge <- rbridge(end = 1, frequency = 500)
  succeeded=all(bridge>=0)||all(bridge<=0)
}
bridge = abs(bridge)
plot(bridge)

wprowadź opis zdjęcia tutaj

RUser4512
źródło
2
t=0
Rzeczywiście, było małe zastrzeżenie;) „A im większa„ częstotliwość ”, którą ustawiłeś, tym mniejsze prawdopodobieństwo znalezienia trajektorii.” ... Jestem tylko w połowie zadowolony z mojej odpowiedzi, ale to jedyna rzecz, o której mogłem pomyśleć o tym, czy musiałbym zacząć od mostu Browna. Szukam (i czekam) na lepszą odpowiedź!
RUser4512,