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 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300
|
\name{Circulant Embedding}
\alias{Circulant}
\alias{RPcirculant}
\alias{Cutoff}
\alias{RPcutoff}
\alias{Intrinsic}
\alias{RPintrinsic}
\title{Circulant Embedding methods}
\description{
Circulant embedding is a fast simulation method for stationary
(possibly anisotropic) Gaussian random
fields on regular grids based on Fourier transformations. It is guaranteed to be an exact method
for covariance functions with finite support, e.g. the spherical
model. The method is admissible for any dimension apart from memory
restrictions. \cr
The simulation is performed on a torus which represents the bended
grid. To remove wrong dependencies occuring at different borders of the
grid which would be close on the torus, the simulation area is
multiplied by a natural number.
There is also a multivariate version of circulant embedding.
Cut-off embedding is a fast simulation method for stationary,
isotropic Gaussian random fields on square lattices based on
the standard \command{\link{RPcirculant}} method,
so that exact simulation is garantueed
for further covariance models, e.g. the \command{\link{RMwhittle}} model.
In fact, the circulant embedding is called with the cutoff
hypermodel, see \command{\link{RMcutoff}}.
\code{Cutoff} halfens the maximum number of
elements models used to define the covariance function of interest
(from 10 to 5).
Here, multiplicative models are not allowed (yet).
\cr
For details see \command{\link{RMcutoff}}.
Intrinsic embedding is a fast simulation method for intrinsically stationary,
isotropic Gaussian random fields on square lattices based on
the standard \command{\link{RPcirculant}} method,
for further \emph{variogram} models, e.g. \command{\link{RMfbm}}.
Note that the simulated random field is always \emph{non-stationary}.
In fact, the circulant embedding is called with the Intrinsic
hypermodel, see \command{\link{RMintrinsic}}.
Here, multiplicative models are not allowed (yet).
\cr
For details see \command{\link{RMintrinsic}}.
}
\usage{
RPcirculant(phi, boxcox, force, mmin, strategy,
maxGB, maxmem, tolIm, tolRe, trials, useprimes, dependent,
approx_step, approx_maxgrid)
RPcutoff(phi, boxcox, force, mmin, strategy,
maxGB, maxmem, tolIm, tolRe, trials, useprimes, dependent,
approx_step, approx_maxgrid, diameter, a)
RPintrinsic(phi, boxcox, force, mmin, strategy,
maxGB, maxmem, tolIm, tolRe, trials, useprimes, dependent,
approx_step, approx_maxgrid, diameter, rawR)
}
\arguments{
\item{phi}{See \command{\link{RPgauss}}.}
% \item{loggauss}{See \command{\link{RPgauss}}}
\item{boxcox}{the one or two parameters of the box cox transformation.
If not given, the globally defined parameters are used.
See \command{\link{RFboxcox}} for details.
}
\item{force}{
Logical. Circulant embedding does not work
if the constructed circulant
matrix has negative eigenvalues. Sometimes it is convenient
to replace all the negative eigenvalues by zero
(\code{force=TRUE}) after \code{trials} number of trials.
Default: \code{FALSE}.
}
\item{mmin}{Scalar or vector, integer if positive.
\code{CE.mmin} determines the initial size of the circulant
matrix. If \code{CE.mmin=0} the minimal starting size is
determined automatically according to the
dimensions of the grid.
If \code{CE.mmin>0} then the absolute starting size is given.
If \code{CE.mmin<0} then the automatically determined
matrix size is multiplied
by \eqn{|\code{CE.mmin}|}; here, \code{CE.mmin} must be smaller
than -1; the
value -1 takes over the minimal starting size.\cr
Note: in any cases, the initial size might be increased according
to \code{CE.useprimes}.\cr
Default: \code{0}.
}
\item{strategy}{Logical. \code{0}: If the circulant
matrix has negative eigenvalues then the
size in each direction is doubled; \cr \code{1}: The size is
enhanced only in
one direction, namely that one where the covariance function has the
largest value at the end point of the grid --- note that
the default value of \code{trials} is probably too small
in that case.
In some cases \code{strategy=0} works better, in other cases
\code{strategy=1}. Just try. Clearly, if the field is
isotropic and a square grid should be simulated, then
\code{strategy=0} is the better choice.
Default: \code{0}.
}
\item{maxGB}{Maximal memory used for the circulant matrix in units of
GB.
If this argument is set then \code{maxmem} is set to MAXINT.
Default: 1.
}
\item{maxmem}{Integer.
maximal number of entries in a row of the circulant matrix.
The total amount of memory needed for the internal calculations
is %about 16 (=2 * sizeof(double))
% times as large as \code{maxmem}
% if \code{\link{RFoptions}}\code{()$Storing=FALSE}, and
32 (= 4 * sizeof(double)) time as large (factor 2 is needed as complex numbers
must be considered for calculating the fft of the covariance matrix;
another factor 2 is needed for storing the simulated result).
%if \code{Storing=TRUE}.
The value of \code{maxmem} must be at least \eqn{2^d} times as large as
the number of points to be simulated. Here, \eqn{d} is the
space dimension. In some cases even much larger.
Note that \code{maxmem} can be used to control the automatic
choice of the simulation algorithm. Namely, in case of huge
circulant matrices, other simulation
methods (TBM) might be faster and might be preferred by the user.
If this argument is set then \code{maxGB} is set to \code{Inf}.
Default: \code{MAXINT}.
}
\item{tolIm}{
If the modulus of the imaginary part is less than
\code{tolIm} then the eigenvalue is always considered as
real (independently of the value of \code{force}).
Default: \code{1E-3}.
}
\item{tolRe}{
Eigenvalues between \code{tolRe} and 0 are always considered as
0 and set 0 (independently of the value of \code{force}).
Default: \code{-1E-7}.
}
\item{trials}{Integer.
A larger circulant matrix is likely to make more eigenvalues
non-negative. If at least one of the thresholds \code{tolRe} and
\code{tolIm} are missed then the matrix size is doubled
according to \code{strategy},
and the matrix is checked again. This procedure is repeated
up to \code{trials - 1} times. If there are still negative
eigenvalues, the simulation method fails if \code{force=FALSE}.
Default: \code{3}.
}
\item{useprimes}{Logical. If \code{FALSE}
the columns of the circulant matrix
have length \eqn{2^k} for some \eqn{k}. Otherwise the algorithm
tries to find a nicely factorizable number close to the size of the
given matrix.
Default: \code{TRUE}.
}
\item{dependent}{Logical. If \code{FALSE}
then independent random fields are created. If \code{TRUE}
then at least 4 non-overlapping rectangles are taken out of the
the expanded grid defined by the circulant matrix.
These simulations are dependent.
See \link{RFoptionsAdvanced} for an example.
See \code{trials} for some
more information on the circulant matrix.
Default: \code{FALSE}.
}
\item{approx_step}{Real value.
It gives the grid size of the approximating grid in case
circulant embedding is used although the points do not lie on a grid.
If \code{NA} then \code{approx_step} is chosen such that
\code{approx_maxgrid} is nearly reached.
Default: \code{NA}.
}
\item{approx_maxgrid}{
It defaults to \code{maxmem}.
}
\item{diameter}{See \command{\link{RMcutoff}} or \command{\link{RMintrinsic}}.}
\item{a}{See \command{\link{RMcutoff}}.}
\item{rawR}{See \command{\link{RMintrinsic}}.}
}
\details{
Here, the algorithms by Dietrich and Newsam are implemented.
}
\value{
An object of class \code{\link[=RMmodel-class]{RMmodel}}.
}
\references{
Circulant Embedding
\itemize{
\item Chan, G. and Wood, A.T.A. (1997)
An Algorithm for Simulating Stationary Gaussian Random Fields.
\emph{Journal of the Royal Statistical Society: Series C (Applied
Statistics)}, \bold{46}, 171--181.
\item Dietrich, C. R. and G. N. Newsam (1993)
A fast and exact method for multidimensional gaussian
stochastic simulations.
\emph{Water Resour. Res.} \bold{29(8)}, 2861--2869.
\item Dietrich, C. R. and G. N. Newsam (1996)
A fast and exact method for multidimensional Gaussian stochastic
simulations: Extension to realizations conditioned on direct and
indirect measurements
\emph{Water Resour. Res.} \bold{32(6)}, 1643--1652.
\item Dietrich, C. R. and Newsam, G. N. (1997)
Fast and Exact Simulation of Stationary Gaussian Processes through
Circulant Embedding of the Covariance Matrix.
\emph{SIAM J. Sci. Comput.} \bold{18}, 1088--1107.
\item Wood, A. T. A. and Chan, G. (1994)
Simulation of Stationary Gaussian Processes in \eqn{[0, 1]^d}.
\emph{Journal of Computational and Graphical Statistics}
\bold{3}, 409--432.
}
Cutoff and Intrinsic
\itemize{
\item Gneiting, T., Sevecikova, H, Percival, D.B., Schlather M.,
Jiang Y. (2006) Fast and Exact Simulation of Large {G}aussian
Lattice Systems in {$R^2$}: Exploring the Limits.
\emph{J. Comput. Graph. Stat.} \bold{15}, 483--501.
\item Stein, M.L. (2002) Fast and exact simulation of fractional
Brownian surfaces. \emph{J. Comput. Graph. Statist.} \bold{11}, 587--599
}
}
\me
\seealso{
\link{Gaussian}, \link{RP}
}
\examples{\dontshow{StartExample()}
RFoptions(seed=0) ## *ANY* simulation will have the random seed 0; set
## RFoptions(seed=NA) to make them all random again
model <- RMstable(s=1, alpha=1.8)
x <- seq(-3,3,0.1)
\dontshow{if (!interactive()){print("x values changed"); x <- -1:1}}
z <- RFsimulate(model=RPcirculant(model), x=x, y=x, n=1)
plot(z)
model <- RMexp(var=10, s=2)
z <- RFsimulate(model=RPcirculant(model), 1:10)
plot(z)
model <- RMfbm(Aniso=diag(c(1,2)), alpha=1.5)
z <- RFsimulate(model, x=1:10, y=1:10)
plot(z)
\dontshow{FinalizeExample()}}
\keyword{methods}
|