Можно ли добавлять на график вариограммы линии/текст/небольшие отметки/персонализированные оси?

Я пытаюсь добавить текст, линии, небольшие деления и/или персонализированные оси на вариограмму с помощью функции plot.variogram. Когда я пытаюсь добавить любой из них, используя традиционные вызовы функций (т. е. text("Text Here")), он возвращает ошибку plot.new has not been called yet, хотя окно построения вариограммы открыто.

Вот мой код:

#v is sample variogram, vmf is fitted model
plot(v, model=vmf, xlim=c(0, 65), ylim=c(0,25), xlab="Distance between Point Pairs (km)",
ylab="Semivariance ((C/km) )", cex.xlab=6.5, cex.ylab=6.5, cex.xaxis=2.5, cex.main=5.5)

#Add a 2 to the y label that is in 10 pt. font so it looks like it is (C/km)^2
par(ps=10, cex=1, cex.main=1)
text(-2, 16, labels=2, srt=90)

#Add lines showing the desired point pair distance and semivariance for the problem
par(new=TRUE, lines(c(53,53),c(0,15),col='red'))
par(new=TRUE, lines(c(0, 53),c(15,15),col='red'))

#Add axis minor tick marks in increments of 5
axis(side=1, at=c(0, 5, 15, 25, 35, 45, 55, 65), labels = NA, tck=-0.01, pos=0) 
axis(side=2, at=c(0, 2.5, 7.5, 12.5, 17.5, 22.5, 25),labels = NA, tck=-0.01, pos=0)

Я попытался «обмануть» R, позвонив:

plot(c(0,65), c(0,25))

а затем запустите код выше. Это позволяет работать традиционным функциям, но, к сожалению, они находятся не в соответствующих местах (т. е. x=5 не находится в позиции 5 по оси x).

Любые рекомендации по лучшим способам «обмануть» R для правильного построения графика? Любые функции, автоматически добавляющие текст, оси и т. д. к графикам вариограмм?

Пожалуйста, дайте мне знать, если есть что-то еще, что вы хотели бы знать.

Спасибо!


r plot kriging personalization gstat
person Jared Smith    schedule 02.04.2014    source источник

Ответы (2)


arrow_upward
1
arrow_downward

Я не уверен, как вы получили свою вариограмму, но используя часть кода и информацию из эта ссылка, мне удалось получить то, что может вам помочь.

Используя пакет geoR, функцию variog вы можете манипулировать plot как обычно.

> sampleV <- 
    read.table(header = TRUE, text = "Station   Av8top      Lat       Lon
  1       60 7.225806 34.13583 -117.9236
  2       69 5.899194 34.17611 -118.3153
  3       72 4.052885 33.82361 -118.1875
  4       74 7.181452 34.19944 -118.5347
  5       75 6.076613 34.06694 -117.7514
  6       84 3.157258 33.92917 -118.2097
  7       85 5.201613 34.01500 -118.0597
  8       87 4.717742 34.06722 -118.2264
  9       88 6.532258 34.08333 -118.1069
  10      89 7.540323 34.38750 -118.5347", row.names = 1)

> library(geoR)
> sampleVMF <- variog(coords = sampleV[,3:4], data = sampleV[,2], 
                      breaks = seq(0, 1.5, length = 11))
> plot(sampleVMF, axes = FALSE,  
       xlab="Distance between Point Pairs (km)",
       ylab="Semivariance ((C/km) )")
> axis(1, at = sampleVMF$u)
> axis(2, at = sampleVMF$v)
> box()
> text(median(sampleVMF$u), median(sampleVMF$v), "Hello world!")
> lines(sampleVMF$u, sampleVMF$v)

введите здесь описание изображения

person Rich Scriven    schedule 02.04.2014
comment
Спасибо за помощь! Ключевым моментом здесь является использование пакета geoR вместо пакета gstat. Команда variog в пакете geoR может быть использована для создания вариограммы, затем команда variofit может использоваться для подгонки вариограммы. - person Jared Smith; 03.04.2014
comment
Мне также интересно узнать, возможно ли добавление этих функций с помощью пакета gstat. Не проблема преобразовать код для использования geoR, просто из любопытства. - person Jared Smith; 03.04.2014

arrow_upward
4
arrow_downward

Ответ Ричарда Скривена хорошо работает при использовании пакета geoR.

При использовании пакета gstat все графики необходимо изменить с помощью команд решетчатой ​​графики (пакет lattice и пакет latticeExtra). Причина ошибки plot.new has not been called yet заключается в использовании базовой графики на решетчатой ​​графике.

Пример модификации сюжета:

Создайте график, используя list, чтобы изменить все входные параметры. scales изменяет оси x и y, и вы можете добавлять деления (at) без labels (см. список осей y). key добавляет легенду.

p1 = plot(v.GF, model=vmf.GF, lwd=2, col.line="black", pch=16, cex=1.5, ylim=c(0,25), xlim=c(0,150), 
      main=list("Gaussian Semivariogram Model for Geothermal Gradient", cex=2),
      xlab=list("Distance between Point Pairs (km)", cex=2), 
      ylab=list(expression("Geothermal Gradient Semivariance (°C/km)"^2), cex=2), 
      scales=list(x=list(at=c(0,25,50,75,100,125,150), labels=c(0,25,50,75,100,125,150)), y=list(at=c(0,2.5,5,7.5,10,12.5,15,17.5,20,22.5,25), labels=c(0,"",5,"",10,"",15,"",20,"",25)), cex=2), 
      key=list(text=list(lab=c("Gaussian Model","Sill","Maximum Interpolation Distance")), space="top", lines=list(col=c("black","black","red"), lwd=2, lty=c(1,2,1)), columns=3, cex=1.5))

Затем начертите и измените область построения, используя trellis.focus и llines для добавления линий и ltext для добавления текста:

trellis.focus("panel",1,1)
plot(p1)
trellis.focus("panel",1,1)
llines(x=c(50,50), y=c(0,14.5), col="red", lwd=2, lty=1)
llines(x=c(0,50), y=c(14.5,14.5), col="red", lwd=2, lty=1)
llines(x=c(0,150), y=c(vmf.GF$psill[2]+vmf.GF$psill[1],vmf.GF$psill[2]+vmf.GF$psill[1]), col="black", lty=2, lwd=2)
ltext(x=12,y=5,"Nugget ~6.0", cex=1.5)
trellis.unfocus()

Вот получившийся график: введите здесь описание изображения

person Jared Smith    schedule 14.07.2015