Essa é uma revisão anterior do documento!
Script by PJ on 16/11/2007
<code R> ## ## Importação dos dados de arquivos ## para o R e para um banco Terralib/MySQL ## - teores de elementos medidos em água e localizações das amostras ## - poligonos com municípios do Paraná ##
## carregando os pacotes que serão usados (aRT) require(aRT)
##Conectar com o SGBD MYSQL (local) con=openConn() con
##Apagar banco já existente if(any(showDbs(con)=="geomedicina")) deleteDb(con, "geomedicina", force=T)
##Criando um novo banco db=createDb(con, "geomedicina") db
## string definindo a projeção proj="+proj=latlong +ellps=GRS67 +towgs84=-66.8700,-4.3700,-38.5200"
## ## Layer 1: municípios ## ## Criar no BD layer para polígonos dos municípios l.pr ← createLayer(db, "Parana", proj=proj) l.pr
##Adicionar ao BD o contorno do Paraná disponível em shapefile (fonte: IBGE) addShape(l.pr, tab="Parana", file="parana_pol.shp", id="CODIGO", length=10)
## verificando status do banco e layer db l.pr
## visualizando polygonos plot(l.pr)
## importando polygonos do banco para o R munic ← getPolygons(l.pr) class(munic)
## Q: qual a diferença em usar os dois comandos abaixo? aRTplot(munic) plot(munic)
## importando atributos dos polygonos do BD para o R municTab ← openTable(l.pr, "Parana") municTab municDt ← getData(municTab) dim(municDt) names(municDt) head(municDt)
## centroids e áreas de cada município center ← getOperation(l.pr, "centroid") center ## Q: como converter para UTM (ou obter um UTM) ???
## Q: area deveria usar UTM para o cálculo? ## Q: área calculada com LatLong é correta? area ← getMetric(l.pr, "area") area
## ## Lendo atributos dos poligonos (munícipios) de outra fonte (DATASUS) ## (serão necessárias operações para compatibilizar com dados do IBGE)
figado ← read.csv("cancerfigadopop2005.csv") head(figado)
names(figado) names(municDt)
## checando e corrigindo nomes dos municípios:
## um município trocou o nome….. ## fix (temporário até obter mapa atualizado do IBGE) ## para arquivo do IBGE (Vila Alta == Alto Paraíso)
figado$municipio ← as.character(figado$municipio) figado[figado$municipio == "Alto Paraíso", "municipio"] ← "Vila Alta"
all.equal(municDt$NOME,figado$municipio) ##estão em ordem diferente… all.equal(sort(municDt$NOME),sort(figado$municipio)) ## ainda 1 incongruencia: foo1 ← which(sort(municDt$NOME) != sort(figado$municipio)) cbind(sort(municDt$NOME),sort(figado$municipio))[foo1,] foo2 ← cbind(sort(municDt$NOME),sort(figado$municipio))[foo1,] foo2 municDt[municDt$NOME == foo2[1], "NOME"] ← foo2[2] ## nomes agora OK! all.equal(sort(municDt$NOME),sort(figado$municipio))
## reordenando segundo a ordem do IBGE figado ← figado[match(municDt$NOME, figado$municipio),] all.equal(municDt$NOME,figado$municipio)
## checando e corrigindo CODIGOS dos municípios: ## código do IBGE tem 1 campo a mais que o do DATASUS figado$id ← as.character(figado$id)
cbind(figado$id, municDt$CODIGO) figado$id ← paste(figado$id,substr(municDt$CODIGO, 7, 7),sep="") all.equal(municDt$CODIGO,figado$id)
## Q: suponha agora que eu tenha varios atributos no R, provenientes de outra ## fonte e associados aos polígonos como acima. ## Qual a melhor maneira de adicionar à tabela já existente (Paraná)? ## (tem que ser coluna a coluna com updatecolumns ou melhor abrir outra tabela ? ## No exemplo abaixo quero adicionar o data frame "figado" ao banco
## tentei pegar is ID's da tabela que veio do banco mas não deu… #> getID(municDt) #Erro em getID(municDt) : Not implemented yet #> getID(municTab) #Erro em getID(municTab) : Not implemented yet
## 1. adicionando à tabela já existente: ## (o que tentei a seguir nao funcionou, tem que ser coluna a coluna…???)
## 1.1 tentando adcionar todas as colunas de uma só vez:
## CODIGO foi usado como ID para munic no addShape() portanto uando aqui tb para figado all.equal(municDt$CODIGO,figado$id) municTab updateColumns(municTab, figado) ## nao funcionou mas pelo vignette aRTtable deveria… ## Erro em if (poskey != 1) { : argumento tem comprimento zero
## 1.2 adicionando as colunas uma a uma… (arrrgh) names(figado) for(cc in names(figado)[3:4]){
tmp <- data.frame(id=figado$id, data = figado[,cc]) names(tmp)[2] <- cc print(head(tmp)) createColumn(municTab, cc, type="int") updateColumns(municTab, tmp)
} ## mesmo erro que acima: ##Erro em if (poskey != 1) { : argumento tem comprimento zero
## 2. adicionando como uma nova tabela
## 2.1 usando creatTable() figadoTab ← createTable(l.pr, "CancerFigado") figadoTab l.pr updateColumns(figadoTab, figado) figadoTab
## usando importTable() figadoTab1 ← importTable(l.pr, "FigadoCancer", id="id", figado)
l.pr
## Q: Pedro. qual a diferença entre usar createTable() e importTable() ? ## A segunda simplesmente encapsula os comandos da primeira em uma unica chamada?
## ## Layer 2: dados geoquímicos ##
## Q: ## Pedro, aqui estou tentando usar 3 formas diferentes de colocar ## o dado de pontos no banco, mas as baseadas im importSpData() ## estao dando problemas
## ## Forma 1: (não funcionou, veja msg de erro abaixo) ## teores ← read.table("aguafix.csv", head=T, sep=";") head(teores) summary(teores) plot(munic) points(teores[,c("LONGITUDE", "LATITUDE")], pch=20, cex=0.3, col="blue")
## convertendo teores para SpatialPointsDataFrame names(teores) coordinates(teores) ← c("LONGITUDE","LATITUDE") teores$ID ← as.character(1:nrow(teores)) ## Q: como atribuir projeção ao Layer criado com importSpData()? O seguinte é válido? attr(teores, "proj4string") ← CRS(proj) l1.geoq ← importSpData(db, teores, "GQ1") ## MSG DE ERRO: #Erro em .aRTcall(object, "cppUpdateColumns", colnames = colnames(data), : # Could not update the table l1.geoq
## ## Forma 2: (não funcionou, veja msg de erro abaixo) ## teores ← read.table("aguafix.csv", head=T, sep=";") teores$ID ← as.character(1:nrow(teores)) #ll ← SpatialPoints(teores[,c("LONGITUDE", "LATITUDE")], proj=CRS(proj)) lldf ← SpatialPointsDataFrame(coords=teores[,c("LONGITUDE", "LATITUDE")],
data=teores[,-(4:5)],
proj = CRS(proj))
l12.geoq ← importSpData(db, lldf, "GQ2") #Erro em .aRTcall(object, "cppUpdateColumns", colnames = colnames(data), : # Could not update the table
## ## Forma 3: (Funcionou) ## ## alternativa a importSpData() seria fazer passo a passo teores ← read.table("aguafix.csv", head=T, sep=";") coordinates(teores) ← c("LONGITUDE","LATITUDE") teores$ID ← as.character(1:nrow(teores)) attr(teores, "proj4string") ← CRS(proj)
l3.geoq ← createLayer(db, "GQ3", proj=proj) l3.geoq addPoints(l3.geoq, teores) l3.geoq tb3.geoq ← importTable(l3.geoq, "Elementos", id="ID", teores) l3.geoq tb3.geoq
## ## Forma 4: (Funcionou, mas…) ## ## alternativa a importSpData() seria fazer passo a passo teores ← read.table("aguafix.csv", head=T, sep=";") teores$ID ← as.character(1:nrow(teores)) lldf ← SpatialPointsDataFrame(coords=teores[,c("LONGITUDE", "LATITUDE")],
data=teores[,-(4:5)],
proj = CRS(proj))
l4.geoq ← createLayer(db, "GQ4", proj=proj) l4.geoq addPoints(l4.geoq, lldf) l4.geoq ## ele nao aceitou usar o nome da tabela "Elementos" já usado acima ## mesmo etando em outro layer !!!!!! tb4.geoq ← importTable(l4.geoq, "Elem", id="ID", lldf) l4.geoq tb4.geoq </code R>
