Skip to content

Instantly share code, notes, and snippets.

/morantest.py Secret

Created July 22, 2014 09:22
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
Star You must be signed in to star a gist
Save anonymous/1d59878309acf4d1cb3a to your computer and use it in GitHub Desktop.
# -*- coding: utf-8 -*-
"""
/***************************************************************************
MoranTest
A QGIS plugin
Moran's I Statistic Test
-------------------
begin : 2014-07-22
copyright : (C) 2014 by test
email : test
***************************************************************************/
/***************************************************************************
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/
"""
# Import the PyQt and QGIS libraries
from PyQt4.QtCore import *
from PyQt4.QtGui import *
from qgis.core import *
# Initialize Qt resources from file resources.py
import resources_rc
# Import the code for the dialog
from morantestdialog import MoranTestDialog
import os.path
import numpy as np
from pysal import W, Moran
import matplotlib.pyplot as plt
class MoranTest:
def __init__(self, iface):
# Save reference to the QGIS interface
self.iface = iface
# initialize plugin directory
self.plugin_dir = os.path.dirname(__file__)
# initialize locale
locale = QSettings().value("locale/userLocale")[0:2]
localePath = os.path.join(self.plugin_dir, 'i18n', 'morantest_{}.qm'.format(locale))
if os.path.exists(localePath):
self.translator = QTranslator()
self.translator.load(localePath)
if qVersion() > '4.3.3':
QCoreApplication.installTranslator(self.translator)
# Create the dialog (after translation) and keep reference
self.dlg = MoranTestDialog()
def initGui(self):
# Create action that will start plugin configuration
self.action = QAction(
QIcon(":/plugins/morantest/icon.png"),
u"Run Test", self.iface.mainWindow())
# connect the action to the run method
self.action.triggered.connect(self.run)
# Add toolbar button and menu item
self.iface.addToolBarIcon(self.action)
self.iface.addPluginToMenu(u"&Moran's I Statistic", self.action)
def unload(self):
# Remove the plugin menu item and icon
self.iface.removePluginMenu(u"&Moran's I Statistic", self.action)
self.iface.removeToolBarIcon(self.action)
# run method that performs all the real work
def run(self):
# 레이어 리스트 채우기
canvas = self.iface.mapCanvas()
layers = canvas.layers()
self.dlg.layer_comboBox.clear()
layerInfo = {}
for layer in layers:
layerName = layer.name()
layerInfo[layerName] = layer
self.dlg.layer_comboBox.addItems(layerInfo.keys())
# show the dialog
self.dlg.show()
# Run the dialog event loop
result = self.dlg.exec_()
# See if OK was pressed
if result == 1:
# do something useful (delete the line containing pass and
# substitute with your code)
# 전역변수 설정
FROM_DIST = int(self.dlg.from_lineEdit.text())
TO_DIST = int(self.dlg.to_lineEdit.text())
BY_DIST = int(self.dlg.by_lineEdit.text())
NAME_FIELD = self.dlg.name_lineEdit.text()
VALUE_FIELD = self.dlg.value_lineEdit.text()
CRITICAL_Z = 1.96
##########################
# 레이어에서 정보 추출
# 레이어 선택
#oLayer = self.iface.activeLayer()
#if not oLayer:
# raise UserWarning(u"레이어를 먼저 선택해야 합니다.") # 종료
oLayer = layerInfo[ self.dlg.layer_comboBox.currentText()]
layerName = oLayer.name()
layerType = oLayer.geometryType()
crs = oLayer.crs()
# ID 리스트 확보
oIDs = oLayer.allFeatureIds()
# Progress 생성
progressMessageBar = self.iface.messageBar().createMessage(u"레이어 정보 수집중...")
progress = QProgressBar()
progress.setMaximum(len(oIDs))
progress.setAlignment(Qt.AlignLeft | Qt.AlignVCenter)
progressMessageBar.layout().addWidget(progress)
self.iface.messageBar().pushWidget(progressMessageBar,self. iface.messageBar().INFO)
# centroid,value(y),name 모으기
centroidList = []
dataList = []
nameList = []
for i, oID in enumerate(oIDs):
progress.setValue(i)
iFeature = oLayer.getFeatures(QgsFeatureRequest(oID)).next()
iGeom = iFeature.geometry().centroid()
centroidList.append(iGeom)
data = iFeature[VALUE_FIELD]
dataList.append(data)
name = iFeature[NAME_FIELD]
nameList.append(name)
# 통계 대상 값 수집
#y = np.array(dataList)
#######################
# FROM_DIST에서 TO_DIST까지 BY_DIST 씩 거리 증가하며 Moran's I 구하기
# 전체 몇 번 돌아가야 하는지 계산
totalCnt, mod = divmod((TO_DIST-FROM_DIST), BY_DIST)
totalCnt += 1
progress.setMaximum(totalCnt)
# 거리를 늘려가며 반복하며 Moran's I 계산
miResults = {}
for i, testDist in enumerate(range(FROM_DIST, (TO_DIST+BY_DIST), BY_DIST)):
progress.setValue(i)
# Weight Matrix 계산 위한 정보 수집
neighbors = {}
weights = {}
yList = []
for i, iID, iCent in zip(range(len(dataList)), oIDs, centroidList):
iRowNeighbors = []
iRowWeights = []
for jID, jCent in zip(oIDs, centroidList):
# 동일 지역인 경우 제외
if iID == jID:
continue
# 기준거리 이내인 경우 인접한 것으로 기록
dist = iCent.distance(jCent)
if dist <= testDist:
iRowNeighbors.append(jID)
iRowWeights.append(1)
# iID 지역에 대한 인접 지역 및 가중치 기록
if len(iRowNeighbors) > 0:
neighbors[iID] = iRowNeighbors
weights[iID] = iRowWeights
yList.append(dataList[i])
# 인접지역과 가중치를 기준으로 현재 testDist의 Weight Matrix 계산
w = W(neighbors, weights)
# 현재 testDist의 Moran's I 계산
y = np.array(yList)
mi = Moran(y, w, two_tailed=False)
# 결과 저장
miResults[testDist] = mi
# Progress 제거
self.iface.messageBar().clearWidgets()
###########################
# 거리에 따른 z 값 변화 그래프
# 그래프를 위한 값 모으기
distList = miResults.keys()
distList.sort()
zList = []
for dist in distList:
zList.append(miResults[dist].z_norm)
# 이미 그래프 창이 있더라도 닫기
plt.close()
# 값 그리기
plt.plot(distList, zList, "b")
# 축 정보
plt.xlabel("Distance = d")
plt.ylabel("Z[d]")
# 유효값 기준선
plt.plot([distList[0], distList[-1]], [CRITICAL_Z, CRITICAL_Z], "r")
plt.text(distList[-1], CRITICAL_Z, " Critical\n Z-Value\n %.3f" % CRITICAL_Z, color="r")
# 제목
plt.title("Z[d]s over a range of search distance")
# 그래프 띄우기
plt.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment