@@ -121,21 +121,35 @@ function Base.rand(rng::AbstractRNG, d::AlphaStable)
121121 bϕ = π/ 2 + β* ϕ
122122 x = 2 / π * (bϕ* tan (ϕ) - β* log (π/ 2 * w* cosϕ/ bϕ))
123123 α == 1 || (x += β * tan (π* α/ 2 ))
124-
124+
125125 return loc + scale* x
126126end
127127
128128
129129"""
130- ...
130+
131+ Generate alpha-sub-Gaussian (aSG) random numbers.
132+
133+ The implementation is based on https://github.com/ahmd-mahm/alpha-SGNm/blob/master/asgn.m
134+
135+ Reference:
136+ A. Mahmood and M. Chitre, "Generating random variates for stable sub-Gaussian processes
137+ with memory", Signal Processing, Volume 131, Pages 271-279, 2017.
138+ (https://doi.org/10.1016/j.sigpro.2016.08.016.)
139+
140+
131141# Arguments
132- - `alpha `: characteristic exponent associated with the aSGN(m) process. This is
142+ - `α `: characteristic exponent associated with the aSGN(m) process. This is
133143a scalar input and should lie within `collect(1.1:0.01:1.98)`.
134144- `R`: covariance matrix of any adjacent `m+1` samples in an aSGN(m) process.
135145The dimension of `R` is equal to `m+1`. It should be a symmetric toeplitz matrix.
136146The maximum acceptable size of `R` is `10x10`
137147- `n`: number of samples required
138- ...
148+
149+ # Examples
150+ ```jldoctest
151+ julia> x = rand(AlphaSubGaussian(n=1000))
152+ ```
139153"""
140154Base. @kwdef struct AlphaSubGaussian{T,M<: AbstractMatrix } <: Distributions.ContinuousUnivariateDistribution
141155 α:: T = 1.50
@@ -184,23 +198,10 @@ function subgausscondprobtabulate(α, x1, x2_ind, invRx1, invR, vjoint, nmin, nm
184198 (1 / sqrt (kappa))* kmarg* vjointR/ vjointR1
185199end
186200
187- """
188- Generate alpha-sub-Gaussian (aSG) random numbers.
189-
190- The implementation is based on https://github.com/ahmd-mahm/alpha-SGNm/blob/master/asgn.m
191-
192- Reference:
193- A. Mahmood and M. Chitre, "Generating random variates for stable sub-Gaussian processes
194- with memory", Signal Processing, Volume 131, Pages 271-279, 2017.
195- (https://doi.org/10.1016/j.sigpro.2016.08.016.)
196201
197- # Examples
198- ```jldoctest
199- julia> x = rand(AlphaSubGaussian(n=1000))
200- ```
201- """
202- function Base. rand (rng:: AbstractRNG , d:: AlphaSubGaussian )
202+ function Random. rand! (rng:: AbstractRNG , d:: AlphaSubGaussian , x:: AbstractArray )
203203 α= d. α; R= d. R; n= d. n
204+ length (x) >= n || throw (ArgumentError (" length of x must be at least n" ))
204205 α ∈ 1.10 : 0.01 : 1.98 || throw (DomainError (α, " α must lie within `1.10:0.01:1.98`" ))
205206 m = size (R, 1 )- 1
206207 funk1 = x -> (2 ^ α)* sin (π* α/ 2 )* gamma ((α+ 2 )/ 2 )* gamma ((α+ x)/ 2 )/ (gamma (x/ 2 )* π* α/ 2 )
0 commit comments