Geoestatística e Scripts em R

Script – Análise estatística no R: Estudo de caso com dados do Lago Walker

A explicação do script a seguir está no artigo Análise estatística no R: Estudo de caso com dados do Lago Walker

GKS1.R

setwd("C:\\GK21\\dados\\walker_data")
dados=read.csv("walker_dat.csv",sep=";",
header=T)
x=dados$Xlocation; y=dados$Ylocation
setwd("C:\\GK21\\figuras")
jpeg("mapaLocal.jpeg",width=5,height=5,
units="in",res=300)
plot(x,y,xlab="LESTE",ylab="NORTE",asp=1)
dev.off()

GKS2.R

setwd("C:\\GK21\\dados\\walker_data")
dados=read.csv("walker_dat.csv",sep=";",
header=T)
V=dados$V
n=length(V)
nc=ceiling(1+3.222*log10(n))
#calculando (nc+1) quebras para nc classes
quebras=c(rep(0,(nc+1)))
for (i in 1:(nc+1)){
quebras[i]=min(V)+(i-1)*(max(V)-min(V))/nc}
#setwd("C:\\GK21\\figuras")
#jpeg("histoV.jpeg",width=5,height=5,
#units="in",res=300)
h=hist(V,breaks=quebras)
#dev.off()
print(h)

GKS3.R

setwd("C:\\GK21\\dados\\walker_data")
dados=read.csv("walker_dat.csv",sep=";",
header=T)
V=dados$V; n=length(V)
quebras=seq(min(V),max(V),length=100)
#setwd("C:\\GK21\\figuras")
#jpeg("acumulaV.jpeg",width=5,height=5,
#units="in",res=300)
h=hist(V,breaks=quebras,plot=F)
freq_acum=100*cumsum(h$counts)/sum(h$counts)
plot(h$mids,freq_acum,xlab="V",
ylab="Frequência acumulada",
main="Distribuição de V")
#dev.off()

GKS4.R

moda_francis=function(x,nc){
  #quando a classe modal for 1 ou nc
  #a moda e o ponto medio
  #altere o nc para melhor definicao
  quebras=c(rep(0,(nc+1)))
  zc=(max(x)-min(x))/nc
  for (i in 1:(nc+1)){
  quebras[i]=min(x)+(i-1)*zc}
  h=hist(x,breaks=quebras,plot=FALSE)
  contagens=h$counts
  maximo<-max(contagens)
  indice<-which.max(contagens)
  largura<-quebras[2]-quebras[1]
  if (indice==1){
  d1=1; d2=1} else {
  if (indice==nc){
  d1=1; d2=1}
  else {
    d1<-contagens[indice]-contagens[indice-1]
    d2<-contagens[indice]-contagens[indice+1]}
  }
  moda<-quebras[indice]+(d1/(d1+d2))*largura
  return(moda)
}
descritivas=function(x){
  library(moments)
  n=length(x)
  media=mean(x)
  x50=median(x)
  vmin=min(x)
  vmax=max(x)
  S2=var(x)
  S=sd(x)
  cv=S/media
  moda=moda_francis(x,12)
  print(moda)
  ass=skewness(x)
  kur=kurtosis(x)
  qua=quantile(x)
  aiq=IQR(x)
  res=c(n,media,x50,moda,S2,S,cv,ass,
  kur,qua,aiq,vmin,vmax)
  return(res)
}
setwd("C:\\GK21\\dados\\walker_data")
dados=read.csv("walker_dat.csv",sep=";",header=T)
V=dados$V[dados$V!=-99.00]; 
#estatisticas para V
res=round(descritivas(V),digits=3)
estats=c(
"No. de dados","Média","Mediana",
"Moda","Variância","Desvio Padrão",
"Coef. de variação",
"Amplitude interquartil",
"Assimetria","Curtose","Valor mínimo",
"Quartil Inferior","Quartil Superior",
"Valor máximo")
valores=c(round(res[1],digits=0),
res[2],res[3],res[4],res[5],
res[6],res[7],res[15],res[8],
res[9],res[16],res[11],res[13],
res[17])
print(cbind(estats,valores))

Arquivos para utilização:

Comments

Deixe um comentário

O seu endereço de e-mail não será publicado. Campos obrigatórios são marcados com *