# Computing the spectrum's posterior distribution Derives the posterior distribution of the spectrum of one or several time series, based on data and prior specifications. ```r bspec(x, ...) ## Default S3 method: bspec(x, priorscale=1, priordf=0, intercept=TRUE, two.sided=FALSE, ...) ``` ## Arguments - `x`: a time series object of the data to be analysed. May be a univariate (`ts` object) or multivariate (`mts` object) time series. - `priorscale`: **either** a `numeric` vector giving the scale parameters of the spectrum's prior distribution; recycled if of length 1. **Or** a `function` of frequency. - `priordf`: **either** a `numeric` vector giving the degrees-of-freedom parameters of the spectrum's prior distribution; recycled if of length 1. **Or** a `function` of frequency. - `intercept`: a `logical` flag indicating whether to include the intercept (zero frequency) term. - `two.sided`: a `logical` flag indicating whether to refer to a one-sided or a two-sided spectrum. In particular affects the interpretation of the prior scale parameters, and sets the default for some methods applied to the resulting `bspec` object via its `two.sided` element. - `...`: currently unused. ## Details Based on the assumptions of a zero mean and a finite spectrum, the posterior distribution of the (discrete) spectrum is derived. The data are modeled using the **Maximum Entropy** (Normal) distribution for the above constraints, and based on the prior information about the spectrum specified in terms of the (conjugate) **scaled inverse $chi-squared$ distribution**. For more details, see the references. ## Returns A list of class `bspec` containing the following elements: - **freq**: a `numeric` vector giving the (Fourier-) frequencies that the spectral parameters correspond to. - **scale**: a `numeric` vector giving the scale parameters of the posterior distributions of the spectral parameters corresponding to the above frequencies. These -internally- always correspond to the **one-sided** spectrum, regardless of the `two.sided` flag (see below). - **df**: a `numeric` vector giving the degrees-of-freedom parameters of the posterior distributions of the spectral parameters corresponding to the above frequencies. - **priorscale**: a `numeric` vector giving the prior scale parameters. - **priordf**: a `numeric` vector giving the prior degrees-of-freedom parameters. - **datassq**: a `numeric` vector giving the sum-of-squares contributed by the data. - **datadf**: a `numeric` vector giving the degrees-of-freedom contributed by the data. - **N**: the sample size of the original time series. - **deltat**: the sampling interval of the original time series. - **deltaf**: the frequency interval of the Fourier-transformed time series. - **start**: the time of the first observation in the original time series. - **call**: an object of class `call` giving the function call that generated the `bspec` object. - **two.sided**: a `logical` flag indicating whether the spectrum is to be interpreted as one-sided or two-sided. ## References Roever, C., Meyer, R., Christensen, N. Modelling coloured residual noise in gravitational-wave signal processing. **Classical and Quantum Gravity**, 28 (1):015010, 2011. tools:::Rd_expr_doi("10.1088/0264-9381/28/1/015010") . See also [arXiv preprint 0804.3853](https://arxiv.org/abs/0804.3853). Roever, C. c("[Degrees-of-freedom\n](https://dcc.ligo.org/LIGO-T1100497/public)", "[ estimation in the Student-t noise model](https://dcc.ligo.org/LIGO-T1100497/public)"). Technical Report LIGO-T1100497, LIGO-Virgo Collaboration, 2011. Roever, C. A Student-t based filter for robust signal detection. **Physical Review D**, 84 (12):122004, 2011. tools:::Rd_expr_doi("10.1103/PhysRevD.84.122004") . See also [arXiv preprint 1109.0442](https://arxiv.org/abs/1109.0442). ## Author(s) Christian Roever, christian.roever@med.uni-goettingen.de ## See Also `expectation`, `quantile.bspec`, `sample.bspec`, `ppsample`, `acf.bspec`, `spectrum` ## Examples ```r # determine spectrum's posterior distribution # (for noninformative prior): lhspec <- bspec(lh) print(lhspec) # show some more details: str(lhspec) # plot 95 percent central intervals and medians: plot(lhspec) # draw and plot a sample from posterior distribution: lines(lhspec$freq, sample(lhspec), type="b", pch=20) ######## # compare the default outputs of "bspec()" and "spectrum()": bspec1 <- bspec(lh) spectrum1 <- spectrum(lh, plot=FALSE) plot(bspec1) lines(spectrum1$freq, spectrum1$spec, col="blue") # (note -among others- the factor 2 difference) # match the outputs: # Need to suppress tapering, padding and de-trending # (see help for "spec.pgram()"): spectrum2 <- spectrum(lh, taper=0, fast=FALSE, detrend=FALSE, plot=FALSE) # Need to drop intercept (zero frequency) term: bspec2 <- bspec(lh, intercept=FALSE) # plot the "spectrum()" output: plot(spectrum2) # draw the "bspec()" scale parameters, adjusted # by the corresponding degrees-of-freedom, # so they correspond to one-sided spectrum: lines(bspec2$freq, bspec2$scale/bspec2$datadf, type="b", col="green", lty="dashed") ######## # handle several time series at once... data(sunspots) # extract three 70-year segments: spots1 <- window(sunspots, 1750, 1819.99) spots2 <- window(sunspots, 1830, 1899.99) spots3 <- window(sunspots, 1910, 1979.99) # align their time scales: tsp(spots3) <- tsp(spots2) <- tsp(spots1) # combine to multivariate time series: spots <- ts.union(spots1, spots2, spots3) # infer spectrum: plot(bspec(spots)) ```

bspec function