library(spatstat) ####################################### ### Vicerozmerne a kotovane procesy ### ####################################### # Priklady # ############ # Dvourozmerny bodovy vzorek: data(amacrine) plot(amacrine) # Koty jsou tridy 'factor' - kategoricke, kvalitativni: amacrine$marks # Muzeme se domnivat, ze 'on' a 'off' se vyskytuji nezavisle na sobe? # Jiny dvourozmerny bodovy vzorek: data(chorley) plot(chorley,chars=c(19,3),cols=c("blue","red")) summary(chorley) # Jsou v ramci zjistenych mist s vyskytem nemoci pripady rakoviny hrtanu rozmisteny nahodne? chorley.extra$plotit() # Sestirozmerny vzorek: data(lansing) plot(lansing) # Kazdy typ zvlast: plot(split(lansing)) # Jsou zavislosti mezi polohami jednotlivych druhu stromu? # Priklad kotovaneho vzorku (kvantitativni koty): data(finpines) plot(unmark(finpines)) # jen polohy bodu - funkce 'unmark' vynecha koty plot(finpines) # i se znazornenymi kotami # Koty jsou spojite; podivame se na ne podrobne: ?finpines marks(finpines) hist(marks(finpines)[,1]) hist(marks(finpines)[,2]) # Ovlivnuje blizkost stromu u sebe jejich vysku? # Operace s kotovanymi vzorky # ############################### # Pomoci funkce 'cut' jsou koty z numerickych hodnot transformovany na faktorove: plot(cut(finpines, "height", breaks=3)) # rozdeli body na 3 typy podle velikosti koty plot(cut(finpines, "diameter", breaks=3)) plot(cut(finpines, breaks=3)) # defaultne bere prvni sloupec kot marks(finpines) # Funkci 'split' lze pouzit pro rozdeleni bodoveho vzorku podle predpisu, # ktery si sami specifikujeme: data(nztrees) plot(nztrees) nndist(nztrees) # vzdalenost od kazdeho bodu k jeho nejblizsimu sousedu cn <- cut(nndist(nztrees),4) # rozdelime na 4 skupiny podle vzdalenosti nejblizsiho souseda plot(split(nztrees,f=cn)) # f musi byt factor, lze vyuzit funkce factor nebo as.factor ### Modely kotovani ### ####################### # Nahodne kotovani (random labelling) - polohy jsou dany, koty vytvoreny nezavisle na polohach. plot(nztrees) koty <- runif(nztrees$n,2,8) # nezavisle rovnomerne koty koty markednztrees <- nztrees %mark% koty # pridame koty k existujicimu bodovemu procesu markednztrees <- setmarks(nztrees,koty) # udela totez plot(markednztrees) # Diskretni koty => ctyrrozmerny bodovy proces: pp <- rMatClust(30,0.1,5) sample(letters[1:4], pp$n, replace=TRUE) factor(sample(letters[1:4], pp$n, replace=TRUE)) mpp_rl <- pp %mark% factor(sample(letters[1:4], pp$n, replace=TRUE)) plot(mpp_rl) # Nahodne slozeni (random superposition) - vicerozmerny bodovy proces vznikne slozenim # nezavislych bodovych procesu. mpp_rs <- superimpose('a' = rMaternII(30,0.05), 'b' = rMaternII(50,0.10), 'c' = rMaternII(70,0.05)) plot(mpp_rs) plot(split(mpp_rs)) summary(mpp_rs) # Nahodne kotovani a nahodne slozeni u Poissonova procesu splyvaji. mppp2 <- rmpoispp(c(3,7),win=square(3),types=c("prvni","druhy")) # dvourozmerny Poissonuv kotovany proces plot(mppp2) summary(mppp2) ### Popisne statistiky ### ########################## # Nejblizsi sousede bodu jednoho typu od bodu druheho typu v dvourozmernem procesu: mppp2prvni <- split(mppp2)$prvni mppp2druhy <- split(mppp2)$druhy N <- nncross(mppp2prvni,mppp2druhy) plot(mppp2,cols=c("blue","red")) arrows(mppp2prvni$x,mppp2prvni$y,mppp2druhy[N$which]$x,mppp2druhy[N$which]$y,length=0.15) # Krizova K-funkce: K_{ij}(r) # lambda_j ... intenzita bodoveho procesu s kotami j # lambda_j K_{ij}(r) ... udava stredni pocet bodu s kotou j # v kouli se stredem v bode s kotou i a polomerem r # Pro i=j (stejne typy bodu) je K_{ii}(r) K-funkce bodoveho procesu s kotami i. # Plati: K_{ij}(r) = K_{ji}(r). # Nahodne kotovani: vsechny krizove K-funkce splyvaji a jsou rovny K-funkci nekotovaneho procesu, # diky spravne normalizaci pomoci lambda_j, resp. lambda. # Nahodne slozeni: krizova K-funkce pro ruzne typy bodu je rovna K-funkci Poissonova procesu # (za predpokladu stacionarity), formalne totiz je # \lambda_j * K_{ij}(r) = E_0^{!i} \Phi_j (b(o,r)). # \Phi_i, \Phi_j jsou ovsem nezavisle, podminena stredni hodnota je obycejnou stredni hodnotou # a je tedy rovna \lambda_j * \pi * r^2. # Odhad krizove K-funkce: plot(Kcross(mpp_rl,i='a',j='b')) plot(Kcross(mpp_rs,i='a',j='b')) plot(Kcross(mppp2,i="prvni",j="druhy")) plot(Kcross(mppp2,i="druhy",j="prvni")) # melo by to vyjit stejne (symetricky)? plot(Kcross(mppp2,i="prvni",j="druhy",correction="translate")) plot(Kcross(mppp2,i="druhy",j="prvni",correction="translate")) # s touto korekci vychazi stejne # Vse v jednom obrazku: plot(alltypes(mpp_rl,"K")) plot(alltypes(mppp2,"K")) plot(alltypes(amacrine,"K")) # Podobne lze vykreslit krizove G-funkce nebo J-funkce. ### Testy nezavislosti ### ########################## # a) nahodna superpozice (K_{ij} = pi*r^2): # Nove realizace simulujeme nahodnym posunutim bodu s kotou "i" a nezavisle na tom nahodnym # posunutim bodu s kotou "j": Ea <- envelope(amacrine, Kcross, nsim=39, i="on", j="off", correction="trans", simulate=expression(rshift(amacrine,radius=0.25))) plot(Ea) # b) nezavisle kotovani (K, K_{ij}, K_{i.} jsou za H0 stejne): Kdif <- function(pp,...,i) { Kidot <- Kdot(pp,...,i=i) K <- Kest(pp,...) dif <- eval.fv(Kidot-K) return(dif) } # Nove realizace simulujeme permutovanim kot pri danych polohach: Eb <- envelope(amacrine, Kdif, nsim=39, i="on", correction="trans", simulate=expression(rlabel(amacrine))) plot(Eb) ### Samostatny (detektivni) ukol: ### ###################################### # Korelacni funkce kot - meri zavislost kot bodu v dane vzdalenosti. # Prozkoumejte napovedu a odhadnete, co dela funkce "markcorr" a jaky vyznam maji jeji hodnoty. plot(markcorr(finpines)) plot(markcorr(mppp2)) M <- markcorr(amacrine, correction="translate", method="density", kernel="epanechnikov", bw=0.02) plot(M)