-
Notifications
You must be signed in to change notification settings - Fork 8
/
image.cov.html
284 lines (230 loc) · 9.33 KB
/
image.cov.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd"><html xmlns="http://www.w3.org/1999/xhtml"><head><title>R: Exponential, Matern and general covariance functions for 2-d...</title>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<link rel="stylesheet" type="text/css" href="R.css" />
</head><body>
<table width="100%" summary="page for image.cov {fields}"><tr><td>image.cov {fields}</td><td style="text-align: right;">R Documentation</td></tr></table>
<h2>
Exponential, Matern and general covariance functions for 2-d
gridded locations.
</h2>
<h3>Description</h3>
<p>Given two sets of locations defined on a 2-d grid efficiently multiplies a
cross covariance with a vector. The intermediate compuations (the setup)
can also be used for fast simulation of the processes on a grid using the circulant
embedding technique.
</p>
<h3>Usage</h3>
<pre>
stationary.image.cov(ind1, ind2, Y, cov.obj = NULL, setup = FALSE,
grid, M=NULL,N=NULL,cov.function="stationary.cov", delta = NULL, cov.args = NULL, ...)
Exp.image.cov(ind1, ind2, Y, cov.obj = NULL, setup = FALSE, grid, ...)
Rad.image.cov(ind1, ind2, Y, cov.obj = NULL, setup = FALSE, grid, ...)
matern.image.cov(ind1, ind2, Y, cov.obj = NULL, setup = FALSE, grid,
M=NULL,N=NULL,theta= 1.0, smoothness=.5)
wendland.image.cov(ind1, ind2, Y, cov.obj = NULL,
setup = FALSE, grid, M = NULL, N = NULL, cov.args=NULL, ...)
</pre>
<h3>Arguments</h3>
<table summary="R argblock">
<tr valign="top"><td><code>ind1</code></td>
<td>
<p>Matrix of indices for first set of locations this is a two column matrix
where each row is the row/column index of the image element. If missing
the default is to use all grid locations.
</p>
</td></tr>
<tr valign="top"><td><code>ind2</code></td>
<td>
<p>Matrix of indices for second set of locations. If missing this is taken to
be ind2. If ind1 is missing ind2 is coerced to be all grid locations.
</p>
</td></tr>
<tr valign="top"><td><code>Y</code></td>
<td>
<p>Vector to multiply by the cross covariance matrix. Y must be the same
locations as those referred to by ind2.
</p>
</td></tr>
<tr valign="top"><td><code>cov.args</code></td>
<td>
<p>Any additional arguments or parameters to the covariance function.</p>
</td></tr>
<tr valign="top"><td><code>cov.obj</code></td>
<td>
<p>A list with the information needed to do the multiplication by
convolutions. This is usually found by using the returned list when
setup=T.
</p>
</td></tr>
<tr valign="top"><td><code>cov.function</code></td>
<td>
<p>Name of the (stationary) covariance function.</p>
</td></tr>
<tr valign="top"><td><code>setup</code></td>
<td>
<p>If true do not do the multiplication but just return the covariance object
required by this function.
</p>
</td></tr>
<tr valign="top"><td><code>delta</code></td>
<td>
<p>A distance that indicates the range of the covariance when it has compact support.
For example this is the theta parameter in the Wendland covariance.</p>
</td></tr>
<tr valign="top"><td><code>grid</code></td>
<td>
<p>A grid list giving the X and Y grids for the image. (See example below.)
This is only required if setup is true.
</p>
</td></tr>
<tr valign="top"><td><code>M</code></td>
<td>
<p>Size of x-grid used to compute multiplication (see notes on image.smooth
for details) by the FFT. If NULL, the default for M is the largest power
of 2
greater than or equal to 2*m where m= length( grid\$x).
This will give an exact
result but smaller values of M will yield an approximate, faster result.
</p>
</td></tr>
<tr valign="top"><td><code>N</code></td>
<td>
<p>Size of y-grid used to compute multiplication by the FFT.</p>
</td></tr>
<tr valign="top"><td><code>theta</code></td>
<td>
<p>Scale parameter for Matern.</p>
</td></tr>
<tr valign="top"><td><code>smoothness</code></td>
<td>
<p>Smoothness parameter for Matern (.5=Exponential)</p>
</td></tr>
<tr valign="top"><td><code>...</code></td>
<td>
<p>Any arguments to pass to the covariance function in setting up the
covariance object. This is only required if setup is TRUE.
For <code>stationary.image.cov</code> one can include <code>V</code> a matrix reflecting a rotation and scaling of
coordinates. See stationary.cov for details. </p>
</td></tr>
</table>
<h3>Details</h3>
<p>This function was provided to do fast computations for large numbers of
spatial locations and supports the conjugate gradient solution in
krig.image. In doing so the observations can be irregular spaced
but their coordinates must be 2-dimensional and be restricted to grid
points.
(The function as.image will take irregular, continuous coordinates and
overlay a grid on them.)
</p>
<p>Returned value: If ind1 and ind2 are matrices where nrow(ind1)=m and
nrow(ind2)=n then the cross covariance matrix, Sigma is an mXn matrix
(i,j) element is the covariance between the grid locations indexed at
ind1[i,] and ind2[j,]. The returned result is Sigma%*%Y. Note that one can
always recover the coordinates themselves by evaluating the grid list at
the indices. e.g. cbind( grid\$x[ ind1[,1]], grid\$y[ind1[,2])) will give
the coordinates associated with ind1. Clearly it is better just to work
with ind1!
</p>
<p>Functional Form: Following the same form as Exp.cov stationary.cov for
irregular locations, the covariance is defined as phi( D.ij) where D.ij is
the Euclidean distance between x1[i,] and x2[j,] but having first been
scaled by theta. Specifically,
</p>
<p>D.ij = sqrt( sum.k (( x1[i,k] - x2[j,k]) /theta[k])**2 ).
</p>
<p>See <code>Matern</code> for the version of phi for the Matern family.
</p>
<p>Note that if theta is a scalar then this defines an isotropic covariance
function.
</p>
<p>Implementation: This function does the multiplication on the full
grid efficiently by a 2-d FFT. The irregular pattern in Y is handled by
padding with zeroes and once that multiplication is done only the
appropriate subset is returned.
</p>
<p>As an example assume that the grid is 100X100 let big.Sigma denote the big
covariance matrix among all grid points ( If the parent grid is 100x100
then big.Sigma is 10K by 10K !) Here are the computing steps:
</p>
<p>temp<- matrix( 0, 100,100)
</p>
<p>temp[ ind2] <- Y
</p>
<p>temp2<- big.Sigma%*% temp
</p>
<p>temp2[ind1]
</p>
<p>Notice how much we pad with zeroes or at the end throw away!
Here the matrix multiplication is effected through convolution/FFT tricks
to avoid creating and multiplying big.Sigma explicitly. It is often faster
to multiply the regular grid and throw away the parts we do not need then
to deal directly with the irregular set of locations.
</p>
<p>Note: In this entire discussion Y is treated as vector. However if
one has complete data then Y can also be interpreted as a image matrix
conformed to correspond to spatial locations. See the last example for
this distinction.
</p>
<h3>Value</h3>
<p>A vector that is the multiplication of the cross covariance matrix with
the vector Y.
</p>
<h3>See Also</h3>
<p>smooth.2d, as.image, krig.image, stationary.cov
</p>
<h3>Examples</h3>
<pre>
# multiply 2-d isotropic exponential with theta=4 by a random vector
junk<- matrix(rnorm(100*100), 100,100)
cov.obj<- stationary.image.cov( setup=TRUE,
grid=list(x=1:100,y=1:100),theta=8)
result<- stationary.image.cov(Y=junk,cov.obj=cov.obj)
image( matrix( result, 100,100)) # NOTE that is also a smoother!
# to do it again, no setup is needed
# e.g.
# junk2<- matrix(rnorm(100**2, 100,100))
# result2<- stationary.image.cov(Y=junk2, cov.obj=cov.obj)
# generate a grid and set of indices based on discretizing the locations
# in the precip dataset
out<-as.image( RMprecip$y, x= RMprecip$x)
ind1<- out$ind
grid<- list( x= out$x, y=out$y)
#
# discretized x locations to use for comparison
xd<- cbind( out$x[ out$ind[,1]], out$y[ out$ind[,2]] )
# setup to create cov.obj for exponential covariance with range= 1.25
cov.obj<- stationary.image.cov( setup=TRUE, grid=grid, theta=1.25)
# multiply covariance matrix by an arbitrary vector
junk<- rnorm(nrow( ind1))
result<- stationary.image.cov( ind1, ind1, Y= junk,cov.obj=cov.obj)
# The brute force way would be
# result<- stationary.cov( xd, xd, theta=1.25, C=junk)
# or
# result<- stationary.cov( xd, xd, theta=1.25) %*% junk
# both of these take much longer
# evaluate the covariance between all grid points and the center grid point
Y<- matrix(0,cov.obj$m, cov.obj$n)
Y[32,32]<- 1
result<- stationary.image.cov( Y= Y,cov.obj=cov.obj)
# covariance surface with respect to the grid point at (32,32)
#
# reshape "vector" as an image
temp<- matrix( result, cov.obj$m,cov.obj$n)
image.plot(cov.obj$grid$x,cov.obj$grid$y, temp)
# or persp( cov.obj$grid$x,cov.obj$grid$y, temp)
# check out the Matern
grid<- list( x= seq(-105,-99,,64), y= seq( 40,45,,64))
cov.obj<- matern.image.cov(
setup=TRUE, grid=grid, theta=.55, smoothness=1.0)
Y<- matrix(0,64,64)
Y[16,16]<- 1
result<- matern.image.cov( Y= Y,cov.obj=cov.obj)
temp<- matrix( result, cov.obj$m,cov.obj$n)
image.plot( cov.obj$grid$x,cov.obj$grid$y, temp)
# Note we have centered at the location (grid$x[16],grid$y[16]) for this case
# using sim.rf to simulate an Matern field
look<- sim.rf( cov.obj)
image.plot( grid$x, grid$y, look)
</pre>
<hr /><div style="text-align: center;">[Package <em>fields</em> version 9.9 <a href="00Index.html">Index</a>]</div>
</body></html>