2014-06-25 23 views
0

我有三个变量(第三个是两个物理值的相对频率)。如何改善与akima的interp?

x <- rnorm(1398) 
y <- rnorm(1398) 
a1 <- rnorm(1398) 
data <- data.frame(x, y, a1) 

第一步是与akima,后者filled.contour。

fld <- with(data, interp(x,y,a1)) 
filled.contour(x=fld$x,y=fld$y,z=fld$z, 
    color.palette=colorRampPalette(c("white", "blue"))) 

这里出现的问题是色标上的范围是0.0008-0.0018。 但是当我检查

> max(a1) 
[1] 0.004291845 
> min(a1) 
[1] 0.0007153076 

如何计算FLD充分代表的数据?

+1

我无法复制此内容。按预期将'range(a1)'作为'-2.931849 3.077989'。从干净的R会话再试一次。 – Spacedman

+0

更具体地说,'set.seed(101);范围(rnorm(1398))给出了-3.177210 3.178489,并且应该总是给出大约(-3,3)的值。否则'a1'正在以不同的方式产生,或'rnorm'被破坏(!) –

+0

是的,它是以其他方式产生的。 – MikiBV

回答

2

interp的文档,关于xo论点:

输出网格的x坐标的矢量。 默认值是在x的范围内均匀间隔的40点 。

所以在fld的X可能不包括相应的最大值或最小值A1值x值。

通过增加xoyo的点数,最大内插值将接近实际最大值。例如,

set.seed(100) 
x <- rnorm(1398) 
y <- rnorm(1398) 
a1 <- rnorm(1398) 
data <- data.frame(x, y, a1) 
fld <- with(data, interp(x,y,a1)) 
fld2 <- with(data, interp(x,y,a1, xo=seq(min(x), max(x), length=1000), yo=seq(min(y), max(y), length=1000))) 

max(a1, na.rm=TRUE) 
#[1] 2.949 
max(fld$z, na.rm=TRUE) 
#[1] 2.481 
max(fld2$z, na.rm=TRUE) 
#[1] 2.902 

此外,如果你想插Z到包括你的最大和最小A1,添加相应的x和y的值XO哟。例如,这是如何让它包含a1的最大值。

max.a1.x <- x[which.max(a1)] 
max.a1.y <- y[which.max(a1)] 
# these have to be sorted, since filled.contour will expect them to be. 
xo <- sort(c(seq(min(x), max(x), length=40), max.a1.x)) 
yo <- sort(c(seq(min(y), max(y), length=40), max.a1.y)) 


fld3 <- with(data, interp(x,y,a1, xo=xo, yo=yo)) 
filled.contour(x=fld3$x,y=fld3$y,z=fld3$z, 
    color.palette=colorRampPalette(c("white", "blue"))) 

max(fld3$z, na.rm=TRUE) 
# [1] 2.949 
+0

虽然它永远不会这么糟糕。你试过了吗? – Spacedman

+0

@Spacedman为什么不呢?是的,我试过了。 –

+0

在原始示例中,a1 < - rnorm(1398)'产生0.0007的“min(a1)”和0.0042的“max(a1)”。永远不会发生。它将成为几个西格玛。 OP做了一些我们看不到的错误。他说这是一个“相对频率”,这让我觉得它不是我们所得到的。 – Spacedman

0

我在pastbin中添加了原始数据。

http://pastebin.com/ydL0fhfD

我认为这可能被理解为data.frame。

关于这个问题,马修建议更多的观点,是的,我试过这个,但仍然只是稍微改善了结果。

> max(fld2$z, na.rm=TRUE) 
[1] 0.002222537 
> max(a1, na.rm=TRUE) 
[1] 0.004291845 
+0

此信息更适合作为您帖子的附录,而不是作为答案。我已更新我的答案,向您展示如何强制interp输出最大/最小a1值。 –

+0

谢谢,现在完美! – MikiBV