In pbiecek/SPAG: Indexes of spatial agglomeration

knitr::opts_chunk$set(echo = TRUE)

Wstęp

Mnogość danych mikro-geograficznych pozwala na stosowanie znacznie szerszego wahlarza metod analizy w celu osiągnięcia znaczących wniosków. Posiadając dane przestrzenne na zadanym przedziale czasu dotyczące pewnego zjawiska, można modelować te wielkości jako pewną realizację stochastycznego procesu punktowego. Takie podejście stosowane było już m. in. w analizie rozwoju lasu, a także stanowi podstawę teoretyczną w dalszej analizie.

TBD

Konstrukcja teoretyczna indeksu

TBD

Opis Kolejnych Indeksów

TBD

Konstrukcja teoretyczna indeksu

Korzystając z postaci indeksu SPAG jako iloczynu trzech członów, można uprościć jego postać:

$$I_{SPAG} = I_{o}\cdot I_{c} \cdot I_{d} = \frac{| \bigcup\limits_{i=1}^{n} P_{i} | }{\sum\limits_{i=1}^{n} |P_i|} \cdot \frac{ \sum\limits_{i=1}^{n} |P_i|}{|P|} \cdot I_{d} = \frac{| \bigcup\limits_{i=1}^{n} P_{i} | }{|P|} \cdot I_{d}$$

interpertując powstawanie i rozwój firm jako pewien proces stochastyczny, a także oznaczając

$$I_{co} = \frac{| \bigcup\limits_{i=1}^{n} P_{i} | }{|P|}$$, otrzymuję, że indeks SPAG jest iloczynem dwóch zmiennych losowych. Załóżmy, że $I_{co}$ ma rozkład beta o gęstości: $$f_{I_{co}}(x) = \frac{x^{a-1}(1-x)^{b-1}}{B(a,b)}$$ zaś $I_{d}$ rozkład gamma o gęstości: $$f_{I_d} (x) = \frac{ \lambda^{\beta} x^{\beta - 1} e^{-\lambda x} }{\Gamma (\beta)}$$

Parametry rozkładu gamma zostaną wyliczone, korzystając z następujących estymatorów:

$$\hat{\beta} = \frac{\mu^2}{\sigma^2}$$ $$\hat{\lambda} = \frac{\mu}{\sigma^2}$$

gdzie $\mu$ i $\sigma^2$ to odpowiednio wartość oczekiwana i wariancja zmiennej losowej $I_d$ o rozkładzie $\Gamma (\beta, \lambda)$.

Parametry rozkładu beta wyznaczane są ze wzorów:

$$\hat{\alpha} = \left( \frac{1-\mu}{\sigma^2}-\frac{1}{\mu} \right) \cdot \mu^2$$

$$\hat{\beta} = \hat{\alpha} \cdot \left(\frac{1}{\mu}-1 \right)$$

Eksperyment

Indeks SPAG można traktować jako zmienną losową, której rozkład zależy od obszaru, na którym jest liczona. W poniższym eksperymencie będzie on liczony dla województwa lubelskiego. W celu otrzymania próbki o wystarczającym rozmiarze dane wejściowe zostały dziesięciokrotnie podzielone na zestawy po 500 firm, na których wywołana została funkcja SPAG.

Zaimplementowana w pakiecie metoda zwraca zarówno wyliczoną wartość $I_{d}$ jak i unię kół będących graficzną interpretacją indeksu SPAG, na podstawie których można wyznaczyć wartości członu $I_{co}$:

setwd("C:/Users/Max/Desktop/MGR_ROZNE/SPAGstarerepo/SPAG/materials/aglomeracja/")
data<-read.csv("geoloc data.csv", header=TRUE, sep=";", dec=".")
data$zatr<-ifelse(data$GR_LPRAC==1, 5, ifelse(data$GR_LPRAC==2, 30, ifelse(data$GR_LPRAC==3,150, ifelse(data$GR_LPRAC==4, 600, 1500))))
TestCompanies <- data[,c(23,24,25,12)]
n <- nrow(TestCompanies)

MapDF <- as.data.frame(MapPoland)
lubelskie <- MapPoland[MapDF$jpt_nazwa_=="lubelskie",]

final <- c()
for(i in 1:10){
  samples <- split(sample(1:n),  rep(1:(n/500+1), each=500))
  ns <- length(samples)

  ee <- sapply(1:ns,function(i){
    t <- TestCompanies[samples[[i]],]
    r <- SPAG(t,lubelskie)
    IDist <- tail(r$IDist,1)
    unionCirc <- attr(r, "unionAreaList")
    return(c(IDist,gArea(unionCirc$Total)))
  })
  ee <- t(ee)
  final <- rbind(final, ee)
}

final[,2] <- final[,2]/gArea(lubelskie)

Do sprawdzenia czy zmienne pochodzą z odpowiednich rozkładów skorzystam z testu Kołmogorova-Smirnova:

load('C:/Users/Max/Desktop/AnalizaRozkladu/final.rda')

varGamma <- var(final[,1])
eGamma <- mean(final[,1])


# dla rozkładu gamma:

wspBeta = eGamma^2/varGamma # shape
wspLambda = eGamma/varGamma # rate

ks.test(final[,1], "pgamma",wspBeta,wspLambda, exact=TRUE)

print("wysokie p-value sugeruje, że rozkład został dobrze dobrany")
hist(final[,1],freq=F, xlab="Distance Index", main="Gamma Distribution")
curve( dgamma(x, shape=wspBeta,rate=wspLambda),0.84, 0.96, add=T, col="blue")

varBeta <- var(final[,2])
eBeta <- mean(final[,2])

alfa = ((1-eBeta)/varBeta-1/eBeta)*eBeta^2
beta = alfa*(1/eBeta-1)

ks.test(final[,2], "pbeta",beta,alfa)
print("niskie p-value sugeruje, że rozkład jest zły")

hist(final[,2],freq=F, ylim=c(0,20), xlab="CO Index", main="Beta Distribution" )
curve( dbeta(x, shape1=alfa, shape2=beta),0.40, 0.6, add=T, col="blue")