text:
# Method¶

This document gives a comprehensive overview over Grond’s methodical background. It describes how the objective function and data weighting are defined, the optimisation algorithm works and how model uncertainties are estimated.

The very core of any optimisation is the evaluation of an objective function or misfit value between observed and predicted data. This is most often based on the difference

${\bf d}_{\mathrm{obs}} - {\bf d}_{\mathrm{synth}},$

but can also be any other comparison, like a correlation measure for instance.

Observed data $${\bf d}_{\mathrm{obs}}$$ are post-processed data (or features) derived from the raw measurements. For example, in the context of seismic source inversion, seismic waveforms are tapered to seismic phases of interest, restituted to displacement and filtered. Predicted data $${\bf d}_{\mathrm{synth}}$$ are then modelled seismograms which are tapered and filtered in the same way as the observed waveforms.

## Forward modelling with pre-calculated Green’s functions¶

The forward modelling of raw synthetic data $${\bf d}_{\mathrm{raw, synth}}$$ for earthquake source models requires the calculation of the Green’s function (GF) between all source points and receiver positions involved, based on a medium (velocity) model. In the general earthquake source problem, the positions of the sources change during the optimisation because the misfit is calculated for many different source-receiver configurations. The calculation of the GFs for each specific source-receiver pair is computationally costly and would be a significant contribution to the total computational duration of an optimisation. Therefore, Grond uses pre-calculated GFs, stored in a database called Pyrocko GF store. Such GF Stores can be created with the Fomosto GF management tool of Pyrocko.

Different applications need different types of GF stores. For the purpose of forward modelling in Grond, we have to distinguish between

1. GFs for dynamic seismic waveforms and
2. GFs for static near-field displacement.

Ready to use GF stores can be found in our online repository. Global GFs for several standard earth models, like e.g. the global 1D PREM model are available, as well as regional distance GFs for many profiles of the CRUST 2.0 earth model database. Custom GF stores can be created using the Fomosto tool and an appropriate choice of the numerical method to calculate the GFs (Fomosto back end).

### GFs for seismic waveforms¶

For regional data analyses the QSEIS method for layered media by Wang et al. (1999) is appropriate. For global forward models the QSSP method also by Wang et al. (2017) is more suited.

### GFs for static near-field displacements (measured by using GNSS or InSAR)¶

For the calculation of purely static coseismic surface displacements the use of the PSGRN/PSCMP method by Wang et al. (2006) is suggested for fast forward modelling.

For more details on GF stores, see the Pyrocko documentation.

## Objective function design¶

The objective function (or misfit function) gives a scalar misfit value how well the source model fits the observed data. A smaller misfit value is better than a large one. The source model that results in the smallest values of the objective function is the global minimum of the misfit function optimum model.

The objective function defines what a model fit is and how good or poor models are scaled with respect to others. Furthermore, the objective function has rules how different observed data sets are handled, which Lp-norm is applied and how data errors are considered in optimisations.

Figure 1: Overview of Grond’s objective function design. Each optimisation target (waveform, satellite and campaign GNSS) handles weights similarly and bootstraps differently. Details on how each target and weight vector is formed is described in the section below.

### Misfit calculation and objective function¶

The core of an optimisation is the data-point-wise calculation of the difference between observed and predicted data:

$|{\bf d}_{\mathrm{obs}} - {\bf d}_{\mathrm{synth}}|.$

Grond supports different seismological observations and a combination of those, thus $${\bf d}_{\mathrm{obs}}$$ and $${\bf d}_{\mathrm{synth}}$$ can be:

• Seismic waveforms
• in time domain
• in spectral domain
• in logarithmic spectral domain
• trace’s spectral ratios
• Static surface displacements
• from unwrapped InSAR images
• from pixel offsets
• measured by using GNSS sensors

The misfit is based on the configurable $$L^p$$-norm with $$p \,\, \epsilon \,\, [1, 2, 3, ...]$$:

(1)$\lVert e \rVert_p = \lVert {\bf{d}}_{\mathrm{obs}} - {{\bf d}}_{\mathrm{ synth}} \rVert_p = \ \left(\sum{|{ d}_{\mathrm{obs}, i} - {d}_{\mathrm{ synth}, i}|^p}\right)^{\frac{1}{p}}$

Further the misfit normalisation factor $$norm$$ is associated with each target. This measure will be used to normalise the misfit values for relative weighting:

(2)$\lVert e_{\mathrm{0}} \rVert_p = \lVert {\bf{d}}_{\mathrm{obs}} \rVert_p = \left(\sum{|{d}_{\mathrm{obs},i}|^p} \right)^{\frac{1}{p}}.$

The resulting normalised misfit

(3)$\lVert e_{\mathrm{norm}} \rVert_p = \ \frac{\lVert e \rVert_p}{ \lVert e_{\mathrm{0}} \rVert_p}.$

is a useful measure to evaluate the data fit. Model predictions that manage to explain parts of the observed data holds $$\lVert e_{\mathrm{norm}} \rVert_p <1$$. Furthermore, the data norm $$\lVert e_{\mathrm{0}} \rVert_p$$ is used in the normalisation of data groups.

#### Waveform misfit¶

Waveform data is preprocessed before misfit calculation: Before the misfit is calculated, observed and synthetic data are tapered within a time window and bandpass filtered (see above). The misfit in Grond can further be based on the maximum waveform correlation.

When measuring waveform data’s cross-correlation, the misfit function is based on the maximum correlation $$\mathrm{max}(C)$$ of $${\bf d}_{\mathrm{obs}}$$ and $${\bf d}_{\mathrm{synth}}$$ defined as:

(4)\begin{align*} e_{\mathrm{cc}} = \frac{1}{2} - \frac{1}{2}\, \mathrm{max}(C), \, \ \mathrm{with} \,\,\, e_{\mathrm{0, cc}} = \frac{1}{2} \,\, ,\, \mathrm{such\,\, that} \ e_{\mathrm{norm}} = 1 - \mathrm{max}(C). \end{align*}

#### Satellite misfit¶

The surface deformation data is pre-processed with kite (See example project: Rectangular source plane from InSAR observations) to obtain a subsampled quadtree. The misfit is then calculated for each quadtree tile $$d_{i}$$.

#### GNSS misfit¶

Each GNSS component (north, east and up) is forward modelled and compared with the observed data.

### Target Weighting¶

Grond implements several different kinds of weights:

• $$w_{\mathrm{tba},i}$$ - target balancing (for waveforms and GNSS campaign only).
• $$w_{\mathrm{noise},i}$$ - noise-based data weights (for waveforms only).
• $$w_{\mathrm{man},i}$$ - user-defined, manual weights of target groups.
• normalisation within data groups, leads to balancing between data groups.

These weights are applied as factors to the misfits, optionally as a product of weight combinations. E.g. for a waveform all data weights combined means:

(5)$w_{\mathrm{comb},i} = w_{\mathrm{tba},i} \cdot w_{\mathrm{noise},i} \ \cdot w_{\mathrm{man},i}$

The misfit and data norm calculations with data weights $$w_{\mathrm{comb},i}$$ change to:

(6)\begin{split}\begin{align*} \lVert e \rVert_x &= \left(\sum{ ({w_{\mathrm{comb},i}} \cdot |{{d}}_{\mathrm{obs},i} - \ {{ d}}_{\mathrm{synth},i}|)^{x}}\right)^{\frac{1}{x}}\\ \lVert e_{\mathrm{0}} \rVert_x &= \left(\sum{ ({w_{\mathrm{comb},i}} \cdot \ |{{d}}_{\mathrm{obs},i} |)^{x}}\right)^{\frac{1}{x}} \end{align*}\end{split}

#### Target balancing weights¶

With these weights waveform targets are balanced with respect to the expected earthquake signal amplitude.

Figure 2: Qualitative sketch how target balancing weight increases with source-receiver distance to balance amplitude inferred by geometrical spreading. Large triangles indicate larger weights and vice versa.

Signal amplitudes in a trace $$|{\bf{d}}_{\mathrm{synth}}|_l$$ depend on the (1) source-receiver distance, (2) on the phase type and (3) the signal processing applied (taper or bandpass). The problem tackled with this particular weight is that large signal amplitudes have higher contributions to the misfit than smaller signals, without providing more information about the source mechanism. From synthetic waveforms of N forward models that have been randomly drawn from the defined model space the mean signal amplitude of the traces is derived. The weight for each trace is then the inverse of these mean signal amplitudes:

(7)${\bf w}_{\mathrm{tba}}=\frac{1}{{\bf a}},\quad \textrm{with trace weights}\,\,a_l= \frac{1}{N} \sum^N_j{|{\bf d}_{\mathrm{synth}}|_{jl}}.$

These balancing weights will enhanced small signals and supress large signals in the objective function. This is described as adaptive station weighting in the PhD thesis by Heimann (2011) (page 23). In Grond they are defined as balancing weights and are received from the TargetBalancingAnalyser module before the optimisation.

#### Data weights based on data error statistics¶

There are direct data weight vectors $$\bf{w}$$ or weight matrices $$\bf{W}$$ based on empirical data error variance estimates. Partly, e.g. for InSAR and GNSS data, these weights are derived from data error correlations expressed in the data error variance-covariance matrix $$\bf{\Sigma}$$:

(8)${\bf w} = \frac{1}{{\bf \sigma}}, \quad \bf{W} = \sqrt{{\bf \Sigma}^{-1}}.$

For a WaveformTargetGroup the data error statistics stem from real data noise before the first phase arrival as described e.g. in Duputel et al. (2012). From the noise traces the inverse of their standard deviation is used. In Grond they are named station_noise_weights and are received from the NoiseAnalyser before the optimisation. In the current grond version, the noise variance and expected target signal (see target balancing above) are multiplied, and not added, which would be more correct. We propose to use the NoiseAnalyser to take out very noisy stations using the options mode='weeding' and a corresponding cutoff factor.

For a SatelliteTargetGroup the data error statistics are pre-calculated by Kite and loaded with the scenes. The estimation of the noise statistics has to be done before Grond by using Kite. In Kite the noise estimation can be done in areas of the displacement map that are not affected by coseismic deformation by using spatial sampling methods and semi-variogram and covariogram formation, described e.g. in Sudhaus and Jonsson (2009).

For a GNSSCampaignTargetGroup the data error statistics are also obtained from the data set. They have to be estimated before and given in the GNSS data YAML-file describing the data set. For details visit the corresponding chapter in the Pyrocko tutorial.

#### Manual data weighting¶

User-defined manual data weights enable an arbitrary weighting of data sets in contrast to balancing of single observations through target balancing and noise-based data weights. No rules apply other than from the user’s rationale. In Grond they are called manual_weight and are given in the configuration file of the targets config.

### Normalisation of data and data groups¶

The normalisation in Grond is applied to data groups that are member of the so called normalisation_family. A normalisation family in Grond can be composed in many ways. However, it is often meaningful to put data of the same kind and with similar weighting schemes into the same normalisation family (see also Fig. 1). This could be P and S waves, or two InSAR data sets. As an explanation some examples are given here:

#### Example 1: Fitting waveforms of P and S waves¶

Let’s say we use the waveform fit in time domain and in spectral domain combined. We then have weighted misfits as in Equation (6) for P waves with $${\bf d}_{\mathrm{obs,Pt}}$$ and $${\bf d}_{\mathrm{synth,Pt}}$$ in time domain and $${\bf d}_{\mathrm{obs,Ps}}$$ and $${\bf d}_{\mathrm{synth,Ps}}$$ in spectral domain. We have also the corresponding weighted misfit norms (see Equation (6)) and the same for S waveforms in time and spectral domain. Let’s also say we are using the $$L^2\,$$-norm.

The waveforms of P and S waves in time domain are of a similar and kind and can, maybe even should, be normalised together. The same may be meaningful for the normalisation of the P and S waves in spectral domain. In our formulations the combined targets form combined vectors, e. g. $${\bf d}_{\mathrm{obs,time}} = ({\bf d}_{\mathrm{obs,Pt}}, {\bf d}_{\mathrm{obs,Ps}})$$.

In Grond we say the time- domain data and the spectral-domain data each belong to a different normalisation_family. Hence, in this example, there are two normalisation families ($$N_{\mathrm{norm\_fams}}=2$$). The global misfit of models that reduce the data variance is constrained to be smaller than 1. The number of normalisation families is accounted for in the global misfit calculation accordingly.

The global misfit for two normalisation families will read:

(9)$\lVert e_{\mathrm{norm,\,global}} \rVert_{2} = \sqrt{ \frac{ \left( \frac{ \lVert e_{\mathrm{time}} \rVert_2}{\lVert \ e_{\mathrm{0,time}} \rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{spectral}} \rVert_2}{\lVert \ e_{\mathrm{0,spectral}} \rVert_2 }\right)^2 }{ \ \left( \frac{ \lVert e_{\mathrm{0,time}} \rVert_2}{\lVert \ e_{\mathrm{0,time}}\rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{0,spectral}} \rVert_2}{\lVert \ e_{\mathrm{0,spectral}}\rVert_2}\right)^2 }} = \ \sqrt{ \frac{ \left( \frac{ \lVert e_{\mathrm{time}} \rVert_2}{ \ \lVert e_{\mathrm{0,time}} \rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{spectral}} \rVert_2}{\lVert \ e_{\mathrm{0,spectral}} \rVert_2 }\right)^2 \ }{ N_{\mathrm{norm\_fams}} }}$

#### Example 2: Fitting waveforms of P waves and static surface displacements¶

Let’s say we use P waveforms in the time domain $${\bf d}_{obs,\mathrm{Pt}}$$. We combine the waveform misfit defined in Equation (6) with the misfit of the maximum waveform correlation defined in Equation (4). Furthermore we use InSAR-measured static surface displacements $${\bf d}_{\mathrm{obs,insar}}$$ and GNSS-measured static surface displacements $${\bf d}_{\mathrm{obs, gnss}}$$. The static surface displacement misfit is defined as in Equation (6). InSAR and GNSS targets should commonly be combined in one normalisation family such that $${\bf d}_{\mathrm{obs,static}} = ({\bf d}_{\mathrm{obs,insar}}, {\bf d}_{\mathrm{obs,gnss}})$$.

The waveform misfits and the correlations, even if the same weights are applied, are measures of a different nature. Also the dynamic waveforms and the static near-field displacements have different relationships to the source parameters. Different normalisation is meaningful. The static surface displacement data themselves should be comparable, even though InSAR and GNSS positing are very different measuring techniques.

In this example we have four target groups that form normalizations families. The global misfit in this example is then:

(10)$\lVert e_{\mathrm{norm,\,global}} \rVert_{2} = \sqrt{ \frac{ \left( \frac{ \lVert e_{\mathrm{time}} \rVert_2}{\lVert \ e_{\mathrm{0,time}} \rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{cc}} \rVert_2}{\lVert \ e_{\mathrm{0,cc}} \rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{static}} \rVert_2}{\lVert \ e_{\mathrm{0,static}} \rVert_2 }\right)^2 }{ \ \left( \frac{ \lVert e_{\mathrm{0,time}} \rVert_2}{\lVert \ e_{\mathrm{0,time}}\rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{0,cc}} \rVert_2}{\lVert \ e_{\mathrm{0,cc}}\rVert_2}\right)^2 +\ \left( \frac{ \lVert e_{\mathrm{0,static}} \rVert_2}{\lVert \ e_{\mathrm{0,static}}\rVert_2}\right)^2 }},$

which is in principle:

(11)$\lVert e_{\mathrm{norm,\,global}} \rVert_{2} = \ \sqrt{ \frac{ \left( \frac{ \lVert e_{\mathrm{time}} \rVert_2}{ \ \lVert e_{\mathrm{0,time}} \rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{cc}} \rVert_2}{ \ \lVert e_{\mathrm{0,cc}} \rVert_2}\right)^2 + \ \left( \frac{ \lVert e_{\mathrm{static}} \rVert_2}{\lVert \ e_{\mathrm{0,static}} \rVert_2 }\right)^2 \ }{ N_{\mathrm{norm\_fams}} }}.$

## The bootstrap method¶

Bootstrapping in Grond (see also Bootstrapping [Wikipedia]) enables to suppress some types of bias in the optimization results. Observations that are affected by other signals or noise often show large misfits. Also insufficient media models for the forward model can result in high misfit values. Already a few high misfit values may pull the optimisation to a biased optimum. With bootstrapping techniques we can better estimate model parameter uncertainties in an efficient way. These include the propagation of the data error, but also the assessment of modelling errors to some extent.

In Grond the bootstrapping is applied in a number of parallel bootstrapping chains where individual bootstrap weights and bootstrap noise is applied to the model misfits. Technically each bootstrap chain carries out its optimization. Find more detail below, at Bayesian Bootstrap Optimisation (BABO). (What is an optimiser?)

In Grond two different bootstrapping methods are implemented:

1. Bayesian and classic bootstrapping through misfit weighting and
2. Residual bootstrapping by adding synthetic noise to the residuals (Fig. 1).

### Bayesian bootstrapping¶

These bootstrap types are based on residual weighting. We divert from the physics-related and noise-related target weights and create numerous additional random weight factors for each target. Virtually equal weights of 1 for each target are redistributed to new random weights, which add up to equal the number of targets. In this way the final misfit values are comparable even without normalisation.

#### Classic weights¶

For a classic bootstrap realisation we draw $$N_{\mathrm{targets}}$$ random integer numbers $${\bf r} \, \in \, [0, N_{\mathrm{targets}}]$$ from a uniform distribution (Fig. 2, left). We then sort these in $$N_{\mathrm{targets}}$$ bins (Fig. 2, right). The frequency in each bin forms the bootstrap target weights.

Figure 3: Formation of classical bootstrap weights. Uniformly random samples (left) and the corresponding histogram (right) with the occurrence frequencies being used as bootstrap weights.

#### Bayesian weights¶

For a Bayesian bootstrap realisation we draw $$N_{\mathrm{targets}}$$ random real numbers $${\bf r} \, \in \, [0, N_{\mathrm{targets}}]$$ from a uniform distribution (Fig. 4, left). We then sort the obtained random values in an ascending order and ensure $$r_1 = 0$$ (Fig. 4, middle). The bootstrap weight now is the difference between two samples (Fig. 4, right):

$w_{\mathrm{bootstr},\,i}=r_{i+1}-r_i$

Figure 4: Formation of Bayesian bootstrap weights. Uniformly random samples (left) are sorted (middle) and the differences of neighbouring points (right) are being used as bootstrap weights.

### Residual bootstrap¶

Residual bootstrapping is a computationally more efficient implementation of the Randomize-then-optimize approach: with empirical estimates of the data error statistics individual realisations of synthetic correlated random noise are systematically added to the data to obtain perturbed optimisations results (Fig. 5).

Earthquake source parameter distributions retrieved with the Randomize-then-optimize method based on the data error variance-covariance matrices have been shown to match the model parameter distributions obtained through Markov Chain Monte Carlo sampling of the model space (Jonsson et al., 2014). In our residual bootstrapping method we add realisations of synthetic correlated random noise to each bootstrapping chain (Fig. 5C and 1). This reduces the calculation of many independent forward models compared to Randomize-then-optimize approach.

To generate random noise we use functions of the Kite module. From the noise estimation region defined in the Kite scenes (Fig. 5A), the noise power spectrum is used directly with a randomised phase spectrum to create new random noise with same spectral characteristics (Fig. 5B). The noise is then subsampled through the same quadtree as defined for the observed data (Fig. 5C).

Figure 5: Residual bootstrapping of InSAR surface displacement data in Grond. (A) From displacements maps we extract noise and (B) synthesise random correlated data noise. (C) This synthetic noise is then subsampled exactly as the observed data. These random realisations are added to the residuals of each bootstrap chain.

## Optimisation¶

Grond’s modular framework is open for different optimisation schemes, the native optimisation schemes is the so-called Bayesian Bootstrap Optimisation (BABO). The Optimiser defines the particular objective function or objective functions and options for them. The optimiser also defines the model space sampling schemes. Multiple objective functions are realized in parallel running optimisation chains - the bootstrap chains (see below).

### Bayesian Bootstrap Optimisation (BABO)¶

Bayesian bootstrap optimisation BABO allows for earthquake source optimisation whilst providing the complete information for a fully Bayesian analysis. BABO is based on Direct Search, where random model parameters are drawn from a defined model space. Those synthetic models are then calculated and compared with the target’s observed data. This needs no assumptions on the topology of the misfit space and is appropriate for highly non-linear problems.

BABO can be configured for a simple Monte-Carlo random direct search. It can also resemble a simulated annealing optimisation approach. Last but not least BABO enables fully probabilistic bootstrapping of the optimisation results. This is realised in parallel with optimisation chains to which bootstrapping weights are applied.

Note

Weights are explained above. The specific weighting is configured with the targets config used and also with the problem. The model space in which the optimisation takes place is defined with the problem. Here described is the sampling and in the context of the multiple objective functions given by the bootstrapping.

#### Sampling scheme and sampling phases¶

Like in any direct search optimisation models are drawn from the model space. From all visited and evaluated models we form and keep a so-called highscore list. The sampling is set up to progressively converge to the low-misfit regions efficiently. However, for multi-modal model parameters distributions an efficient sampling can loose sight of multiple minima with significantly low misfits. In Grond we can use measures to nurse these multiple minima.

##### Highscore list¶

This list contains a defined number of the current best models (lowest misfit). It is continuously updated at runtime. The highscore list length $$L_{\mathrm{hs}}$$ (i.e. number of member models) is problem dependend: $$L_{\mathrm{hs}} = f_{\mathrm{len}} (N_{\mathrm{par}} -1)$$, with $$N_{\mathrm{par}}$$ being the number of model parameters. $$f_{\mathrm{len}}$$ is configurable (chain_length_factor, default is 8).

There are three sampling phases defined, based on which models are drawn from the model space:

Figure 7: Strategic sketch of different optimiser sampling phases.

##### UniformSamplerPhase¶

At the beginning of the optimisation this sampler phase explores the solution space uniformly. A configurable number of models are drawn randomly from the entire model space based on a uniform distribution.

##### InjectionSamplerPhase¶

This starting phase allows to inject pre-defined models at the start of the optimisation. These models could originate from a previous optimisation.

##### DirectedSamplerPhase¶

This sampler is used for the second phase and follows any of starting samplers above: Using existing models of the current highscore models the directed sampler draws a configurable number of new models. Like this convergence to low-misfit regions is enabled. There are quite some noteworthy configurable details to this sampler phase:

###### sampling_distributions¶

New models are drawn from normal distribution. The standard deviations are derived from the highscore models parameter’s standard deviation and scaled by scatter_scale (see below). Optionally, the covariance of model parameters is taken into account by configuring when multivariate_normal is enabled (default is normal distribution). The distribution is centred around

1. mean of the highscore model parameter distributions
2. a random model from the highscore list or
3. an excentricity_compensated draw (see below).
###### scatter_scale¶

This scales search radius around the current highscore models. With a scatter scale of 2 the search for new models has a distribution with twice the standard deviation as estimated from the current highscore list. It is possible to define a beginning scatter scale and an ending scatter scale. This leads to a confining directed search. In other words, the sampling evolves from being more explorative to being more exploitive in the end.

###### starting_point¶

This method tunes to the center value of the sampler distribution: This option, will increase the likelihood to draw a highscore member model off-center to the mean value. The probability of drawing a model from the highscore list is derived from distances the highscore models have to other highscore models in the model parameter space. Eccentricity is therefore compensated, because models with few neighbours at larger distances have an increased likelihood to be drawn.

What’s the use? Convergence is slowed down, yes, but to the benefit of low-misfit region represented by only a few models drawn up to the current point.

Let’s assume there are two separated groups of low-misfit models in our highscore list, with one group forming the 75% majority. In the directed sampler phase the choices of a mean center point for the distribution as well as a random starting point for the sampler distribution would favour new samples in the region of the highscore model majority. Models in the low-misfit region may be dying out in the highscore list due to favour and related sparse sampling. eccentricity compensations can help is these cases and keep models with not significantly higher misfits in the game and in sight.

TODO: correct? too many explanations? Sebastian, here is the perfect place for one of your movies.

### Bootstrap chains¶

A bootstrap chain is a realisation of target bootstrap weights and/or target bootstrap residuals (depending on the targets, Fig. 7A). Therefore each bootstrap chain has a different misfit. With one forward model $$N_{\mathrm{bootstrap}}$$ different global misfits are calculated (Fig. 7B). This approach allows many bootstrap chains leeching the same forward models.

The highscore list member models in each bootstrap chain (Fig. 7B) will differ to some extent and therefore different bootstrap chains may converge to different places within the model space (Fig. 7C, Fig. 8). These differences mark the uncertainty of the models with respect to data errors.

Figure 7: Bootstrap chain graph. (A) Illustration of bootstrap weights, larger triangles indicate higher weight. (B) bootstrap chain highscore lists and (C) their influence on the convergence in the model parameter space due to the individual weighted objective function of each bootstrap chain.

The convergence of model parameters for the models within each bootstrap chain is dependent on the settings of the optimisation, e.g. the setup of parameter bounds, scatter scale settings of the directive sampling phase and other tunable. With very exploitive settings convergence can be forced. However, if the convergence within each bootstrap chain starts to form individual clusters in the model space, further optimisation will not provide significantly better models. In Fig. 8 the area of the highscore models of three bootstrap chains has only little overlap compared to an earlier stage visualised in Fig. 7C.

Figure 8: Drawing new model candidates from the described sampling strategies - the proposal is based on the existing solution space.

## Movies: BABO at work¶

TODO: replace draft text with something meaningful. add figure with toy problem setup.

Toy problem: find best fitting source location in 3D, given noisy 1D distance measures from 10 observers on the horizontal z=0 plane. Projection to vertical cross section is shown. Star is true solution. Lines indicate regions of low misfit.

### Single chain¶

Only upper half-space is searched, problem is unimodal.

### Bimodal, standard deviations from high score models¶

Full space is searched, problem is bimodal due to symmetry


text:
Under the auspices of the Computational Complexity Foundation (CCF)

REPORTS > DETAIL:

### Revision(s):

Revision #3 to TR16-053 | 12th December 2017 23:57

#### Completeness for First-Order Properties on Sparse Structures with Algorithmic Applications

Revision #3
Authors: Jiawei Gao, Russell Impagliazzo, Antonina Kolokolova, Ryan Williams
Accepted on: 12th December 2017 23:57
Keywords:

Abstract:

Properties definable in first-order logic are algorithmically interesting for both theoretical and pragmatic reasons. Many of the most studied algorithmic problems, such as Hitting Set and Orthogonal Vectors, are first-order, and the first-order properties naturally arise as relational database queries. A relatively straightforward algorithm for evaluating a property with k+1 quantifiers takes time $O(m^k)$ and, assuming the Strong Exponential Time Hypothesis (SETH), some such properties require $O(m^{k-\epsilon})$ time for any $\epsilon > 0$. (Here, m represents the size of the input structure, i.e. the number of tuples in all relations.)

We give algorithms for every first-order property that improves this upper bound to $m^k/2^{\Theta(\sqrt{\log n})}$, i.e., an improvement by a factor more than any poly-log, but less than the polynomial required to refute SETH. Moreover, we show that further improvement is equivalent to improving algorithms for sparse instances of the well-studied Orthogonal Vectors problem. Surprisingly, both results are obtained by showing completeness of the Sparse Orthogonal Vectors problem for the class of first-order properties under fine-grained reductions. To obtain improved algorithms, we apply the fast Orthogonal Vectors algorithm of [AWY15,CW16].

While fine-grained reductions (reductions that closely preserve the conjectured complexities of problems) have been used to relate the hardness of disparate specific problems both within P and beyond, this is the first such completeness result for a standard complexity class.

Changes to previous version:

Revision #2 to TR16-053 | 21st February 2017 20:58

#### Completeness for First-Order Properties on Sparse Structures with Algorithmic Applications

Revision #2
Authors: Jiawei Gao, Russell Impagliazzo, Antonina Kolokolova, Ryan Williams
Accepted on: 21st February 2017 20:58
Keywords:

Abstract:

Properties definable in first-order logic are algorithmically interesting for both theoretical and pragmatic reasons. Many of the most studied algorithmic problems, such as Hitting Set and Orthogonal Vectors, are first-order, and the first-order properties naturally arise as relational database queries. A relatively straightforward algorithm for evaluating a property with k+1 quantifiers takes time $O(m^k)$ and, assuming the Strong Exponential Time Hypothesis (SETH), some such properties require $O(m^{k-\epsilon})$ time for any $\epsilon > 0$. (Here, m represents the size of the input structure, i.e. the number of tuples in all relations.)

We give algorithms for every first-order property that improves this upper bound to $m^k/2^{\Theta(\sqrt{\log n})}$, i.e., an improvement by a factor more than any poly-log, but less than the polynomial required to refute SETH. Moreover, we show that further improvement is equivalent to improving algorithms for sparse instances of the well-studied Orthogonal Vectors problem. Surprisingly, both results are obtained by showing completeness of the Sparse Orthogonal Vectors problem for the class of first-order properties under fine-grained reductions. To obtain improved algorithms, we apply the fast Orthogonal Vectors algorithm of [AWY15,CW16].

While fine-grained reductions (reductions that closely preserve the conjectured complexities of problems) have been used to relate the hardness of disparate specific problems both within P and beyond, this is the first such completeness result for a standard complexity class.

Changes to previous version:

(1) Algorithmic applications: Our reduction from first-order properties to OV gives better algorithms for first-order properties.

(2) Extending the result to hypergraphs (relations of arity greater than 2).

(3) The new reduction algorithm is deterministic.

(4) We are glad that Antonina Kolokolova and Ryan Williams became our coauthors.

Revision #1 to TR16-053 | 27th April 2016 03:10

#### Orthogonal Vectors is hard for first-order properties on sparse graphs

Revision #1
Authors: Jiawei Gao, Russell Impagliazzo
Accepted on: 27th April 2016 03:10
Keywords:

Abstract:

Fine-grained reductions, introduced by Vassilevska-Williams and Williams, preserve any improvement in the known algorithms. These have been used very successfully in relating the exact complexities of a wide range of problems, from NP-complete problems like SAT to important quadratic time solvable problems within P such as Edit Distance. However, until now, there have been few equivalences between problems and in particular, no problems that were complete for natural classes under fine-grained reductions. We give the first such completeness results. We consider the class of first-order graph property problems, viewing the input in adjacency list format (aka "sparse graph representation"). For this class, we show that the sparse Orthogonal Vectors problem is complete under randomized fine-grained reductions. In proving completeness for this problem, we also show that this sparse problem is equivalent to the standard Orthogonal Vectors problem when the number of dimensions is polynomially related to the number of vectors. Finally, we also establish a completeness and hardness result for k-Orthogonal Vectors.

Our results imply that the conjecture "not every first-order graph problem has an improved algorithm" is a useful intermediary between SETH and the conjectured hardness of problems such as Edit Distance. It follows that, if Edit Distance has a substantially subquadratic algorithm, then every first order graph problem has an improved algorithm. On the other hand, if first order graph property problems have improved algorithms, this falsifies SETH (and even some weaker versions of SETH) and gives new circuit lower bounds. We hope that this is the beginning of extending fine-grained complexity to include classes of problems as well as individual problems.

Changes to previous version:

Updated the k=2 case in Section 5, "Step 1-2".

### Paper:

TR16-053 | 6th April 2016 00:29

#### Orthogonal Vectors is hard for first-order properties on sparse graphs

TR16-053
Authors: Jiawei Gao, Russell Impagliazzo
Publication: 8th April 2016 14:56

format: plain


text:
# I Direction of rotational frame dragging

Tags:
1. Nov 15, 2017

### Sonderval

I have the following question considering frame dragging:
A test mass starting at rest near a rotating mass or with an initial velocity pointing towards the center of the rotating mass will be deflected in such a way that it begins to move around the mass in the rotational direction. This is normal frame dragging.
Now consider the opposite situation: A test mass starts near the rotating mass but is moving outwards with a velocity in radial direction. In which direction would this mass be deflected?

According to standard formulations of gravitoelectromagnetism like here
http://www.arxiv.org/abs/gr-qc/0311030
the acceleration is proportional to v×B (with B being the analogue of the magnetic field); so if I reverse the velocity, the acceleration should be reversed as well, i.e., if a moving mass is deflected to the right relative to its velocity on falling in, a mass moving out should also be deflected to the right relative to its velocity.

But when I try to see this effect (for example using the program GRorbits
http://stuleja.org/grorbits/), the outwards moving mass is deflected in the co-rotating direction.

So I'm confused which answer is correct and where the mistake in my argument concerning the cross-product lies, if there is one.

I suspect that part of the problem may be that the "outward radial direction" depends on whether I look at the situation from the point of view of a distant Schwarzschild observer or a local frame, but I am not exactly sure how to see this more clearly or make this idea exact.

2. Nov 16, 2017

### Paul Colby

If I read (I'm using the word read rather loosely here, more like a momentary glance) the paper linked to correctly, $B$ should be along $J$, the angular momentum of the mass. In this case the rotation of $v$ is about $J$ which means the change in $v$ would have the opposite sign if the direction is reversed.

3. Nov 17, 2017

### Sonderval

@Paul Arveson
Thanks. Yes, and this means that if a particle flying in one direction is accelerated to the right, a particle flying in the opposite direction is also accelerated to the right.
Actually, I think I found the problem: I looked only at the deflection of the orbit in GRorbits - but I should look at how the deflection changes over time. This shows that the direction of the acceleration is actually opposite to the rotation of the mass as the formula suggests.
The initial deflection of the orbit is due to the fact that the "outward direction" is not the same for an observer close to the mass as for a distant Schwarzschild observer.


text:
## Van der Pol oscillator

 Flow* vanderpol.model

## Model description

The dynamics of the Van der Pol oscillator can be defined by the following ODE

wherein is a constant whose typical value is 1.

## Reachability settings

We consider the initial set which is used in [1].

## Result

As an example, the following figure shows an overapproximation, computed by Flow*, of the reachable set in the time horizon .

## References

[1] M. Althoff, O. Stursberg, M. Buss. Reachability analysis of nonlinear systems with uncertain parameters using conservative linearization. In Proceedings of the 47th IEEE Conference on Decision and Control (CDC’08), pages 4042 – 4048, IEEE, 2008.


text:
# Revision history [back]

The problem seems to be that X is a generator for the finite field on 9 elements, $\text{GF}(9)$, whereas R is just the integers mod 9, $\mathbb{Z}_9$.

Perhaps what you want is

sage: K(X+5)

The problem seems to be that X is a generator for the finite field on 9 elements, $\text{GF}(9)$, whereas R is just polynomial ring over the integers mod 9, $\mathbb{Z}_9$.$\mathbb{Z}_9[x]$.

Perhaps what you want is

sage: K(X+5)

The problem seems to be that X is a generator for the finite field on 9 elements, $\text{GF}(9)$, whereas R is polynomial ring over the integers mod 9, $\mathbb{Z}_9[x]$.Sorry my previous answer was completely wrong.

Perhaps what if you want isuse

sage: K(X+5)
y=R.gen()

and then

sage: R(y+5)
X+5


text:
# A drone

A flying drone aimed the area for an architect. He took off perpendicularly from point C to point D. He was at a height of 300 m above the plane of ABC. The drone from point D pointed at a BDC angle of 43°. Calculate the distance between points C and B in meters.

Result

x =  279.755 m

#### Solution:

$CD = 300 \ m \ \\ \angle BDC = 43 \ ^\circ \ \\ \angle DCB = 90 \ ^\circ \ \\ \ \\ δ = \angle BDC ^\circ \rightarrow rad = \angle BDC ^\circ \cdot \ \dfrac{ \pi }{ 180 } \ = 43 ^\circ \cdot \ \dfrac{ 3.1415926 }{ 180 } \ = 0.75049 \ \\ \ \\ \tan δ = BC / CD \ \\ \tan δ = x / CD \ \\ \ \\ x = CD \cdot \ \tan ( δ ) = 300 \cdot \ \tan ( 0.7505 ) \doteq 279.7545 = 279.755 \ \text{ m }$

Our examples were largely sent or created by pupils and students themselves. Therefore, we would be pleased if you could send us any errors you found, spelling mistakes, or rephasing the example. Thank you!

Leave us a comment of this math problem and its solution (i.e. if it is still somewhat unclear...):

Be the first to comment!

Tips to related online calculators
Do you want to convert length units?

## Next similar math problems:

1. Angles of elevation
From points A and B on level ground, the angles of elevation of the top of a building are 25° and 37° respectively. If |AB| = 57m, calculate, to the nearest meter, the distances of the top of the building from A and B if they are both on the same side of t
2. Depth angle
From a cliff of 150 meters high, we can see the ship at a depth angle of 9° at sea. How far is the ship from the cliff?
3. Depth angles
At the top of the mountain stands a castle, which has a tower 30 meters high. We see the crossroad in the valley from the top of the tower and heel at depth angles of 32° 50 'and 30° 10'. How high is the top of the mountain above the crossroad
4. Ruler
How far from Peter stands 2m hight John? Petr is looking to John over ruler that keeps at arm's distant 60 cm from the eye and on the ruler John measured the height of 15 mm.
5. Maple
Maple peak is visible from a distance 3 m from the trunk from a height of 1.8 m at angle 62°. Determine the height of the maple.
6. Thales
Thales is 1 m from the hole. The eyes are 150 cm above the ground and look into the hole with a diameter of 120 cm as shown. Calculate the depth of the hole.
7. How far
From the top of a lighthouse 145 ft above sea level, the angle of depression of a boat 29°. How far is the boat from the lighthouse?
8. The Eiffel Tower
The top of the Eiffel Tower is seen from a distance of 600 meters at an angle of 30 degrees. Find the tower height.
9. If the
If the tangent of an angle of a right angled triangle is 0.8. Then its longest side is. .. .
10. Aircraft
The plane flies at altitude 6500 m. At the time of first measurement was to see the elevation angle of 21° and second measurement of the elevation angle of 46°. Calculate the distance the plane flew between the two measurements.


text:
Home About the Journal Latest Work Current Issue Archive Special Issues Editorial Board

2021,  3 (3):   248 - 260

Published Date：2021-6-20 DOI: 10.1016/j.vrih.2021.05.004

Abstract

Background
There is a large group of deaf-mutes worldwide, and sign language is a major communication tool in this community. It is necessary for deaf-mutes to be able to communicate with others who are capable of hearing, and hearing people also need to understand sign language, which produces a great demand for sign language tuition. Even though there have already been a large number of books written about sign language, it is inefficient to learn sign language through reading alone, and the same can be said on watching videos. To solve this problem, we developed a smartphone-based interactive Chinese sign language teaching system that facilitates sign language learning.
Methods
The system provides a learner with some learning modes and captures the learner's actions using the front camera of the smartphone. At present, the system provides a vocabulary set with 1000 frequently used words, and the learner can evaluate his/her sign action by subjective or objective comparison. In the mode of word recognition, the users can play any word within the vocabulary, and the system will return the top three retrieved candidates; thus, it can remind the learners what the sign is.
Results
This system provides interactive learning to enable a user to efficiently learn sign language. The system adopts an algorithm based on point cloud recognition to evaluate a user's sign and costs about 700ms of inference time for each sample, which meets the real-time requirements.
Conclusion
This interactive learning system decreases the communication barriers between deaf-mutes and hearing people.

Content

1 Introduction
The World Health Organization (WHO) reported that there are approximately 466 million people worldwide with disabling hearing loss 1 . This is over 5% of the world's population, and 34 million of these people are children. The Second China National Sample Survey on Disability that the number of people with hearing disabilities in China was about 20.04 million, accounting for 1.53% of the total population at that time[1]. In such a wide community, sign language is a major communication tool that expresses specific meanings with movement, position, and orientation of hands, body, and facial expressions. Sign language is a fully-fledged natural language with its own lexicon and grammar.
With the development of society, an increasing number of people attach importance to learning sign language. However, the teaching of sign language suffers from several challenges, such as the number of qualified teachers, lack of interaction during in-class learning, and limited attention to each student[2]. As smartphones are very popular today, this work aims to provide an interactive teaching system to those who are willing to learn sign language.
In this study, we designed and implemented a Chinese sign language teaching system on a smartphone that covers the Chinese sign language vocabulary set with 1000 frequently used words. As shown in Figure 1, the system provides multiple functions for learning: users can learn the description of a sign language word, sketch, and demonstration video in word learning mode; interactive learning through demonstration show-action play and assessment in learning test mode; and text and sign searching in dictionary searching mode. Compared to in-class sign language learning, this system provides an opportunity for personalized one-to-one learning. The developed system can be a valuable educational tool for out-of-school practice. It can also be a learning platform for those who want to learn sign language and for deaf people who want to learn Mandarin.
2 Related work
2.1 Computer vision for assistive technologies
Over the last few decades, there has been a tremendous increase in demand for Assistive Technologies (AT), which are useful for overcoming functional limitations of individuals and improving their quality of life[3]. Lee and Medioni presented a vision-based wearable RGB-D camera indoor navigation system[4]. The system guides the visually impaired user from one location to another without providing map or GPS information in advance. Kayukawa et al. developed an assistive suitcase system, BBeep, to support blind people when walking through crowded environments[5].
Some studies have also focused on social medical service systems. Zeng et al. developed a hand gesture system for Human Computer Interaction (HCI) and applied it to intelligent wheelchair control[6]. Parmar et al. collected skeletal data of patients with cerebral palsy during rehabilitation using the Kinect device and evaluated the quality of their movements[7]. Feng et al. proposed a vision-based fall detection method for monitoring elderly people in a home care environment[8]. Martinel et al. proposed a food recognition system that can automatically identify the type and quantity of food and estimate its calorie content[9].
2.2 Existing sign language teaching systems
At present, computer assistive technology is used in sign language teaching, which greatly improves learning efficiency. Some sign language teaching systems have also been developed.
Based on a comprehensive understanding of the way and basic requirements of sign language teaching, Zhang designed a sign language teaching platform using a Kinect sensor[10]. The system develops new functions to input and recognize sign language, which can provide support for current sign language teaching.
Sagawa and Takeuchi developed a Japanese sign language teaching system to help sign language studies[11]. The system includes sign language recognition, which can recognize the inputted Japanese sign language gesture and translate it into Japanese, and sign language generation, which can translate Japanese into Japanese sign language and display it in the form of 3DCG (three-dimensional computer graphics) animation.
Phan et al. designed a sign language teaching system called My Interactive Coach (MIC) [12]. This system can teach learners Australian signs, allow them to practice, record their practiced performance, and provide different forms of visual feedback on how well they performed the signs.
2.3 Sign language recognition
Sign language recognition has attracted an increasing number of researchers[13,14]. Since 2013, some research teams have conducted a series of studies on sign language recognition of isolated words based on Convolutional Neural Networks (CNN). Pigou et al. used a CNN to extract single-frame gesture features, and finally connected three fully connected layers for classification[15]. Their method achieved an accuracy of 91.7 % for vocabulary set of 20 Italian sign language words. Köpüklü et al. fused motion information into static images, sent the spatiotemporal features after fusion into the subsequent CNN network for recognition[16].
The CNN only encodes features within the single-frame image data. 3D-CNN can simultaneously model the spatio-temporal features of a slip video to capture the motion information from multiple frames and realize a more global sense. The 3D-CNN was first proposed by Ji et al. and was used for behavior recognition[17]. To achieve gesture recognition with large-scale videos, Li et al. proposed a C3D model based on RGB-D data[18,19]. Huang et al. proposed a type of sign language recognition method based on the spatial-temporal attention mechanism of a 3D-CNN network, which integrated spatial attention into the network, focused on the hand features of interest, and finally used the temporal attention mechanism to classify the important sign language movements[20].
Compared to the two networks mentioned above, Rerrent Neural Network (RNN) is a neural network used to process sequential data, and it is better at capturing long-term contextual information. Therefore, in recent years, more and more research has been based on sign language recognition using RNN. Cate et al. recognized 95 categories of isolated words in sign language by using sequential modeling of RNN speech video features[21]. Chai et al. proposed a Two streams Recurrent Neural Networks (2s-RNN) in 2016, which ranked first in the Chalearn gesture recognition challenge in the same year[22]. Bantupalli and Xie proposed an RNN network based on an inception model for sign language recognition[23]. Liao et al. sent segmented hand features and original RGB data into the BLSTM and achieved accurate identification of sign language isolated words[24].
3 System design and implementation
In this study, we designed and implemented a Chinese sign language teaching interactive app on a smartphone. Users can interactively learn sign language and evaluate their sign action expressions through this system. As shown in Figure 2, the main interface of this app contains three functional modules. Users can enter the word-learning interface by clicking the upper word learning function, with three learning modes, and obtain the learning statistics from this part. After learning words, users can evaluate the quality of the learning process by clicking icons in the lower half of the interface. The bottom tab provides auxiliary functions, such as text searching and sign searching. We explain each of these functional modules in the following sections.
Figure 1 shows the total architecture of the proposed system. We adopt multilevel interfaces to enable users to easily start and switch between interfaces. To provide better feedback for the learning process, we adopt the sign language quality assessment method to score the sign language performed by users in time and the sign language recognition method to search for relevant signs in the database to enable understanding and comparison. Relevant algorithms are introduced in Section 3.4.
3.1 Word learning
This system contains a multimedia vocabulary set of the 1000 most frequently used words in "Basic Gestures of Chinese Sign Language"[25]. To better organize the learning schedule for such a large vocabulary, the app provides several learning modes, including sequential learning, unit learning, and random learning. In the sequential learning, users can learn sign language words step-by-step, and the app will record the learning history and start from the last word reviewed in the most recent session. In unit learning, we divide the 1000 sign words into 31 units, with approximately 30 words per unit and group them into different categories, such as appellation and occupation. Learners can choose the units that they prefer to learn systematically. In contrast, to weaken the connection between words and avoid the occurrence of confusion, we also provide a random learning mode, which provides the learner with a random order of words.
Furthermore, we add a 'Daily learning' function to the app's home page. Learners can set the number of words they learn each day so that they can learn in their personalized schedule. In addition, the proportion of new words could be set to balance the study and review.
Once a learning mode is selected, the user enters the word-demonstration interface. As shown in Figure 3a, the learning interface displays the demonstration video, sketch, and textual description from top to bottom. Learner can obtain a comprehensive understanding of this sign word, and the app also provides pause and multiple playback speeds to enable the user to capture video details at his/her preferred speed.
In particular, we have a function similar to listening and repeating in spoken language learning: Play with me, as shown in Figure 3b. Learners can enter this function by clicking on a video at the interface. The application will first play a demo video and then wait for 3s before recording. Learners can follow the demo video and practice the same word. After recording, the application will rate learners' performance through the sign language quality assessment algorithm, which provides timely feedback to learners to achieve a better learning experience.
3.2 Learning testing
In addition to learning modules, we also provide examination modules so that learners can evaluate their learning results. Since sign language is a method of communication, we will test learners' learning results from two aspects: recognition of sign words and quality of sign performance. Figure 4 shows these two modules.
The identity module mainly examines the learner's understanding of sign language. The app will play a sign video and provide multiple choices, and learners need to choose the correct one from three options. If learners make the wrong choice, the correct answer will be displayed, and learners will understand what they have misunderstood. After that, learners can click the '>' button to generate the next question.
In the Sign module, the app uses an assessment algorithm to evaluate the learner's actions. Learners need to make the corresponding sign action, according to the hint from the app. After that, the app then gives the learner a score based on his or her performance, and learners can keep practicing to achieve higher scores.
3.3 Dictionary searching
The dictionary searching mode is a function of auxiliary teaching, which is convenient for users to quickly find the information they need. This module includes two functions: text searching and sign searching.
As Figure 5a shows, text searching is the retrieval from text to sign words. When a user forgets the action of a sign word, the user can input the word he/she wants to query in the search box at the top of the interface, and this app will return the words using fuzzy word matching.
Sign searching can recognize the inputted signed gesture and translate it into Mandarin. When the user sees a strange sign word, he/she can use this function to query it. As Figure 5b shows, the user clicks the BEGIN button to start recording the video and clicks the button again after he/she finishes the sign action. Subsequently, the recognition algorithm in the system processes the recorded video. After a very short wait, the app will provide the top three candidates.
3.4 Recognition and evaluation algorithm
We use the front camera of the smartphone to capture images, which are used to recognize and evaluate learners' sign actions. First, the point cloud sequence is extracted from the depth video and sent to the improved FlickerNet[26] to extract the spatio-temporal features. Then, we proposed a time sequence partitioning model to enhance core gesture information by block prediction. Finally, all local features
$h$
are concatenated to form video feature
$H$
. The framework of the proposed approach is illustrated in Figure 6. The operators used in the network are supported by mainstream mobile deep learning, which can be deployed easily and quickly.
3.4.1 Point cloud formation
Considering the smartphone's limited computational power and the inability to handle complex network structures such as 3D convolution, we use a point cloud instead of an RGB image as the input. Each pixel in the depth image represents the distance from the camera and the confidence level. When using sign language for communication, the hands mainly appear in front of the body. Therefore, we can separate the hand area by setting a threshold. Subsequently,
$N$
points are sampled (128 points in default) from the hand area through random sampling or farthest point sampling as the initial point cloud input. An example of this is shown in Figure 7.
3.4.2 Network structure
The input of the network is a point cloud of size
$T × N × d$
, where
$T$
is the number of frames,
$N$
is the number of points sampled per frame, and
$d$
is the number of channels per point. We use only the
as our point input for simplicity and efficiency
$( d = 4 )$
. We used a modified version of FlickerNet to extract features. For a sign language video, we obtain a feature map of
$T × N ' × C$
, where
$N '$
is the number of points after down-sampling, and
$C$
is the output dimension of the points.
RNN is not suitable for running on smartphones because of its relatively complex network structure. Each Chinese sign language word usually contains two or three atomic actions. Therefore, we can divide a sign language video into several small segments and make predictions separately. However, it is not realistic to accurately divide each atomic action. Following the idea of PCB[27], we divide the features along the time dimensions to obtain
$p$
time blocks, as shown in Figure 8. For each part, we adopt an average pooling operation to obtain
$p$
local features
$g$
, and use
$1 × 1$
convolution to reduce the dimension of
$g$
, which extracts the new local features
$h$
. Next, each feature
$h$
is input into a fully connected layer to predict the class probability of the input.
During training, we optimize the network by minimizing the sum of the cross-entropy losses of
$p$
classifiers. In the testing phase,
$p$
local feature vectors
$h$
can be connected to form the feature vector
$H$
of the entire video. For sign language recognition, top-3 classes are selected as candidates.
To enable word learning and measure learning progress, we use the algorithm of quality evaluation. Our evaluation criterion is whether the sign can be recognized correctly. The learner signs the word, and if the algorithm can correctly recognize it, the learner's score is higher. Conversely, if the algorithm recognizes it as another word, the learner's score is lower. Specifically, we obtain the probability of recognition for each word using the previous network. These probabilities are normalized to 0-1 using the softmax function. The probability of selecting the target category was magnified 100 times to obtain the final percentage score.
3.4.3 Deployment model
As this is an interactive application with the critical requirement of low latency, it is extremely important that this app can run on an edge device in real time. At present, mainstream deep learning frameworks support the deployment of mobile terminals. PyTorch released PyTorch Mobile at the end of 2019, which supports an end-to-end workflow from Python to deployment on iOS and Android. Users can seamlessly transition from a training model to a deployment model without the need for other tools.
Therefore, we use PyTorch to write this model and train on a desktop, and then convert it to a TorchScript file by calling the " torch.jit.trace" method. Finally, the file in the "assert" folder is placed on Android. In the Android project, PyTorch Android is added as a gradle dependency in build.gradle to the project. Then, the app can read the model in an assert folder.
For the smartphone, we use the HUAWEI Mate20 Pro because it is equipped with a 3D depth perception camera. It uses structural light to capture 3D data in the range of 70cm to obtain depth images. We also tested the application on Huawei Mate30 Pro.
3.5 Interactive design
To be concise, most of the buttons in the system are represented by ICONs, which avoids the use of words and reduces the barriers between languages. The user can click on the ICONs to select and jump to the selected functions. In addition, the application adds sliding controls, such as left-slide and right-slide, to switch words (Figure 9a), and up-slide to return.
When a user records videos, the optimal distance from the smartphone is approximately 50cm. In this case, it is not convenient to reach out and click the button on the screen. Therefore, we added a space gesture control function that can help the user operate the app at a distance.
The app captures the depth image from the front camera of HUAWEI Mate20 Pro and crops the hand area by presetting the threshold, then draws a "red dot" on the screen to indicate the position of the hand area in the entire depth image, as shown in Figure 9b. The learner can control the movement of the "red dot" by moving his/her hand. To satisfy daily habits, the left and right movements are mirrored. If the "red dot" stays in the ICON area for more than 1s, the app will determine that the learner has clicked it. In this way, learners can switch words during the learning process, start and end video recordings, etc. They can also choose whether to turn on this function in the setting interface.
4 Experiments
4.1 Evaluation of recognition
We perform comprehensive studies on the Chinese Sign Language Dataset, collected from 30 invited deaf-mutes. For the common most 100 words in Chinese sign language, each person performed them twice, and 5892 valid videos were obtained. We select the data of 27 people (a total of 5297 videos) as the training set, and the other (a total of 595 videos) as the test set. We train the models from scratch for 200 epochs with a mini-batch size of 32. Adam, with a momentum of 0.9, is used with a learning rate of 10-4, which was divided by 10 at epoch 100.
Experiments on the number of frames. The effective number of frames for each video is between 20 and 70 frames. To solve the problem of video length inconsistency, we uniformly sample T frames from a depth video. As the number of frames increased, the recognition accuracy also increased. However, when increasing from 16 to 32 frames, the performance improvement is not as obvious. The experimental results are presented in Figure 10a.
Experiments on the number of points. We need to extract
$N$
points from the depth image as a hand point cloud. The greater the number of point clouds, the more accurate the description of the hand shape, and the better the recognition performance. The experimental results are presented in Figure 10b. As 128 points and 256 points have similar performance, we use 128 points for subsequent experiments.
Experiments on timing partitioning number. In the temporal block model, the feature vector feature of the video is divided into
$p$
local feature vectors according to the time dimension. We conduct experiments to determine the value of
$p$
. To ensure that it can be divided evenly, the value of
$p$
is a factor of 16 frames of video length. When
$p = 1$
, it is equivalent to no blocking.
$p = 16$
indicates that each frame is processed separately. The final experimental results are presented in Figure 10c. When
$p = 4$
, the best effect is achieved.
To verify the generalizability of our model, we conduct experiments on the SHREC'17 Dataset[28], which is a public dynamic hand gesture dataset. Table 1 details our performance comparison with several state-of-the-art methods. Because the SHREC'17 dataset provides skeletal data of the hand and fingers, most previous works use skeleton sequences as inputs that provide relatively accurate hand pose structures and joint trajectories. The work[29] demonstrated that clouds based methods can achieve excellent performance even if accurate hand poses are difficult to obtain in real-world cases. In contrast, the network structure proposed in this paper is more concise and easier to deploy to a mobile terminal.
Performance comparison (%) on SHREC17 dataset
Method Modality 14 gestures 28 gestures
Key frames[28] Depth sequence 82.9 71.9
SoCJ+HoHD+HoWR[30] skeleton 88.2 81.9
Res-TCN[31] skeleton 91.1 87.3
STA-Res-TCN[31] skeleton 93.6 90.7
ST-GCN[32] skeleton 92.7 87.7
DG-STA[33] skeleton 94.4 90.7
PointLSTM-middle[29] point clouds 95.9 94.7
baseline point clouds 91.6 89.2
Ours point clouds 96.3 95.0
Recognition algorithms play an important role in the implemented system. Through the algorithm, we can not only recognize sign language and provide more accurate candidate results for the word recognition module, but also score more accurately in quality evaluation to provide more effective feedback to users.
4.2 User study
To evaluate the usability of this system, we invited 10 students, who had little sign language learning experience, to take the test. For comparison, we chose two free sign language apps[34,35] in the Android application market to join the test. We first introduced the usage and main functions of each application to students. Each student was allowed to use the software freely for ten minutes, after which they rated the app according to the following five aspects: comprehensibility of learning content, learning feedback, learning interest, interactive experience, and system stability.
Figure 11 shows a statistical graph of the user experience in comparison with the other two experimental apps. The data of each indicator in the figure are the average values after the statistics. Each item is scored from 0 to 5. The higher the project score, the better the user experience is. As can be seen from the figure, compared with the previous two systems, the system proposed in this paper provides feedback to users and promotes more efficient learning. Meanwhile, as our system features an interactive experience, the users felt more engaged.
At the same time, students generally believe that the system can effectively improve the learning and the make the process of learning sign language more interesting, unlike the less engaging learning method of reading books. However, there is a disadvantage in that the operation is a little tedious, and the whole process sometimes requires more time than referring to books.
5 Discussion
We developed a smartphone-based Chinese sign language teaching system with different learning modes. The teaching system integrates recent sign language technologies, which achieve both accurate and efficient results for sign language recognition and quality evaluation. The developed system helps users learn sign language through a smartphone anytime and anywhere, which reduces the communication barriers between the deaf and those not hard of hearing. This system can also be applied to similar teaching scenarios, such as the gestures of traffic police.
This system still has some aspects that can be improved. We will continue to improve functions and interactions so that the system can be used to learn not only isolated words but also continuous sign language. In addition, we use recognition probability as the basis for sign language assessment. Other assessment methods could be explored in future studies. We hope that this app will be applied to public welfare to improve the study and life of the disabled.

Reference

1.

Gong C, Bin C, Lei Z. Study on the prevention and strategy of disability in China. Procedia-Social and Behavioral Sciences, 2010, 2(5): 6906–6913 DOI:10.1016/j.sbspro.2010.05.041

2.

Lu W. Research and suggestions on sign language teaching in Chinese universities. Journal of Changchun University, 2014

3.

Leo M, Medioni G, Trivedi M, Kanade T, Farinella G M. Computer vision for assistive technologies. Computer Vision and Image Understanding, 2017, 154: 1–15 DOI:10.1016/j.cviu.2016.09.001

4.

Lee Y H, Medioni G. RGB-D camera based wearable navigation system for the visually impaired. Computer Vision and Image Understanding, 2016, 149: 3–20 DOI:10.1016/j.cviu.2016.03.019

5.

Kayukawa S, Higuchi K, Guerreiro J, Morishima S, Sato Y, Kitani K, Asakawa C. BBeep: a sonic collision avoidance system for blind travellers and nearby pedestrians. In: Proceedings of the 2019 CHI Conference on Human Factors in Computing Systems. Glasgow Scotland Uk, New York, NY, USA, ACM, 2019, 1–12 DOI:10.1145/3290605.3300282

6.

Zeng J H, Sun Y R, Wang F. A natural hand gesture system for intelligent human-computer interaction and medical assistance. In: 2012 Third Global Congress on Intelligent Systems. Wuhan, China, IEEE, 2012, 382–385 DOI:10.1109/gcis.2012.60

7.

Parmar P, Morris B T. Measuring the quality of exercises. In: 2016 38th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC). Orlando, FL, USA, IEEE, 2016, 2241–2244 DOI:10.1109/embc.2016.7591175

8.

Feng W G, Liu R, Zhu M. Fall detection for elderly person care in a vision-based home surveillance environment using a monocular camera. Signal, Image and Video Processing, 2014, 8(6): 1129–1138 DOI:10.1007/s11760-014-0645-4

9.

Martinel N, Piciarelli C, Micheloni C. A supervised extreme learning committee for food recognition. Computer Vision and Image Understanding, 2016, 148: 67–86 DOI:10.1016/j.cviu.2016.01.012

10.

Zhang S. Design research in Kinect based sign-language teaching system. East China Normal University. 2014

11.

Juan-les-Pins, France, York New, Press ACM, 2002, 137–145 DOI:10.1145/641007.641035

12.

Phan H D, Ellis K, Dorin A. MIC, an interactive sign language teaching system. In: Proceedings of the 30th Australian Conference on Computer-Human Interaction. Melbourne Australia, New York, NY, USA, ACM, 2018, 544–547 DOI:10.1145/3292147.3292237

13.

IEEE transactions on pattern analysis and machine intelligence-cover. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2004, 26(5): 0_1–673 DOI:10.1109/tpami.2004.1273910

14.

Jiang X W, Satapathy S C, Yang L X, Wang S H, Zhang Y D. A survey on artificial intelligence in Chinese sign language recognition. Arabian Journal for Science and Engineering, 2020, 45(12): 9859–9894 DOI:10.1007/s13369-020-04758-2

15.

Pigou L, Dieleman S, Kindermans P J, Schrauwen B. Sign language recognition using convolutional neural networks. European Conference on Computer Vision. Springer, Cham, 2014, 572–578 DOI:10.1007/978-3-319-16178-5_40

16.

Köpüklü O, Köse N, Rigoll G. Motion fused frames: data level fusion strategy for hand gesture recognition. In: 2018 IEEE/CVF Conference on Computer Vision and Pattern Recognition Workshops (CVPRW). Salt Lake City, UT, USA, IEEE, 2018, 2184–21848 DOI:10.1109/cvprw.2018.00284

17.

IEEE transactions on pattern analysis and machine intelligence information for authors. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2015, 37(6): C3 DOI:10.1109/tpami.2015.2427753

18.

Li Y N, Miao Q G, Tian K, Fan Y Y, Xu X, Li R, Song J F. Large-scale gesture recognition with a fusion of RGB-D data based on the C3D model. In: 2016 23rd International Conference on Pattern Recognition (ICPR). Cancun, Mexico, IEEE, 2016, 25–30 DOI:10.1109/icpr.2016.7899602

19.

Li Y N, Miao Q G, Tian K, Fan Y Y, Xu X, Li R, Song J F. Large-scale gesture recognition with a fusion of RGB-D data based on saliency theory and C3D model. IEEE Transactions on Circuits and Systems for Video Technology, 2018, 28(10): 2956–2964 DOI:10.1109/tcsvt.2017.2749509

20.

Huang J, Zhou W G, Li H Q, Li W P. Attention-based 3D-CNNs for large-vocabulary sign language recognition. IEEE Transactions on Circuits and Systems for Video Technology, 2019, 29(9): 2822–2832 DOI:10.1109/tcsvt.2018.2870740

21.

Cate H, Dalvi F, Hussain Z. Sign language recognition using temporal classification. 2017

22.

Chai X J, Liu Z P, Yin F, Liu Z, Chen X L. Two streams recurrent neural networks for large-scale continuous gesture recognition. In: 2016 23rd International Conference on Pattern Recognition (ICPR). Cancun, Mexico, IEEE, 2016, 31–36 DOI:10.1109/icpr.2016.7899603

23.

Bantupalli K, Xie Y. American sign language recognition using deep learning and computer vision. In: 2018 IEEE International Conference on Big Data (Big Data). Seattle, WA, USA, IEEE, 2018, 4896–4899 DOI:10.1109/bigdata.2018.8622141

24.

Liao Y Q, Xiong P W, Min W D, Min W Q, Lu J H. Dynamic sign language recognition based on video sequence with BLSTM-3D residual networks. IEEE Access, 2019, 7: 38044–38054 DOI:10.1109/access.2019.2904749

25.

General Administration of Quality Supervision, Inspection and Quarantine of the People's Republic of China. Basic Chinese Sign Language Gestures. China Standards Press, 2009

26.

Min Y, Chai X, Zhao L, Chen X. FlickerNet: Adaptive 3D gesture recognition from sparse point clouds. British Machine Vision Conference. 2019, 105

27.

Sun Y, Zheng L, Yang Y, Tian Q, Wang S. Beyond part models: Person retrieval with refined part pooling (and a strong convolutional baseline). Proceedings of the European Conference on Computer Vision. 2018, 480–496 DOI:10.1007/978-3-030-01225-0_30

28.

SunY F, Zheng L, Yang Y, Tian Q, Wang S J. Beyond part models: person retrieval with refined part pooling (and a strong convolutional baseline). Computer Vision-ECCV 2018, 2018, 480–496 DOI:10.1007/978-3-030-01225-0_30

29.

de Smedt Q, Wannous H, -P Vandeborre J, Guerry J, le Saux B, Filliat D. 3D hand gesture recognition using a depth and skeletal dataset: SHREC'17 track. 3Dor '17: Proceedings of the Workshop on 3D Object Retrieval. 2017, 33–38 DOI:10.2312/3dor.20171049

30.

Min Y C, Zhang Y X, Chai X J, Chen X L. An efficient PointLSTM for point clouds based gesture recognition. In: 2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR). Seattle, WA, USA, IEEE, 2020, 5760–5769 DOI:10.1109/cvpr42600.2020.00580

31.

De Smedt Q, Wannous H, Vandeborre J P. Skeleton-based dynamic hand gesture recognition. In: 2016 IEEE Conference on Computer Vision and Pattern Recognition Workshops (CVPRW). Las Vegas, NV, USA, IEEE, 2016, 1206–1214 DOI:10.1109/cvprw.2016.153

32.

Hou J X, Wang G J, Chen X H, Xue J H, Zhu R, Yang H Z. Spatial-temporal attention res-TCN for skeleton-based dynamic hand gesture recognition. Computer Vision-ECCV 2018 Workshops, 2019 DOI:10.1007/978-3-030-11024-6_18

33.

Peng W, Hong X P, Chen H Y, Zhao G Y. Learning graph convolutional network for skeleton-based human action recognition by neural searching. Proceedings of the AAAI Conference on Artificial Intelligence, 2020, 34(3): 2669–2676 DOI:10.1609/aaai.v34i03.5652

34.

Chen Y X, Zhao L, Peng X, Yuan J B, Metaxas D N. Construct dynamic graphs for hand gesture recognition via spatial-temporal attention. 2019

35.

Luo Y. Love Sign Language (Version 1.3.5), 2020

36.

Changsha Keenbow Information Technology Co., Ltd. 2020, Keenbow Sign Language (Version 1.4)


text:
# NMST539 | Lab Session 4

## Wishart Distribution

### LS 2017 | Tuesday 21/03/2017

###### Rmd file (UTF8 encoding)

A user-friendly interface (one of many): RStudio.

Manuals and introduction into R (in Czech or English):

• Bína, V., Komárek, A. a Komárková, L.: Jak na jazyk R. (PDF súbor)
• Komárek, A.: Základy práce s R. (PDF súbor)
• Kulich, M.: Velmi stručný úvod do R. (PDF súbor)
• De Vries, A. a Meys, J.: R for Dummies. (ISBN-13: 978-1119055808)

#### Wishart distribution

The Wishart distribution is, as an analogy to the $$\chi^2$$ distribution and the variance inference in the univariate case, a very useful tool for the analysis of the variance-covariance matrix of same random sample $$X_{1}, \dots, X_{n}$$ for $$n \in \mathbb{N}$$ where each $$X_{i}$$ in the sample is a $$p$$-dimensional random vector ($$p$$ different covariates recorded on each subject). Thus, it is probability distributions defined over symmetric and nonnegative-definite random matrices.

In general, the Wishart distribution is a multi-variate generalization of the $$\chi^2$$ distribution: for instance, for a normaly distributed random sample $$\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n}$$ drawn from some multivariate distribution $$N_{p}(\boldsymbol{0}, \Sigma)$$, with the zero mean vector and $$\Sigma$$ being some variance-covariance matrix we have that $\mathbb{X}^{\top}\mathbb{X} \sim W_{p}(\Sigma, n),$ for $$\mathbb{X} = (\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n})^{\top}$$. Similarly, for a random sample $$X_{1}, \dots, X_{n}$$ drawn from $$N(0,1)$$ we have that $\boldsymbol{X}^{\top}\boldsymbol{X} \sim \chi_{n}^2,$ where $$\boldsymbol{X} = (X_{1}, \dots, X_{n})^\top$$. Thus, for a sample of size $$n \in \mathbb{N}$$ drawn from $$N(0,1)$$ the corresponding distribution of $$\mathbb{X}^\top\mathbb{X}$$ is $$W_{1}(1, n) \equiv \chi_{n}^2$$.

The corresponding density function of the Wishart distribution takes the form $f(\mathcal{X}) = \frac{1}{2^{np/2} |\Sigma|^{n/2} \Gamma_{p}(\frac{n}{2})} \cdot |\mathcal{X}|^{\frac{n - p - 1}{2}} e^{-(1/2) tr(\Sigma^{-1}\mathcal{X})},$ for $$\mathcal{X}$$ being a $$p\times p$$ random matrix and $$\Gamma_{p}(\cdot)$$ is a multivariate generalization of the Gamma function $$\Gamma(\cdot)$$ In the R software there are various options (packages) on how to use and apply the Wishart distribution.

#### To Do by Yourself (Theoretical and Practical)

• Let the multivariate random sample $$\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n}$$ comes from some general multivariate distribution $$N_{p}(\boldsymbol{\mu}, \Sigma)$$, with some mean vector $$\boldsymbol{\mu} \in \mathbb{R}^p$$ and the variance-covariance positive definite matrix $$\Sigma$$. Show, that the random matrix defined as $$n\mathbb{X}^\top \mathcal{H}\mathbb{X}$$ follows the Wishart distribution $$W_{p}(\Sigma, n - 1)$$, for $$\mathcal{H}$$ being the idenpotent centering matrix $$\mathcal{H} = \mathbb{I}_{n} - \frac{1}{n}\boldsymbol{1}_{n}\boldsymbol{1}_{n}^\top$$.
• use the help session in R to find out more about some standard commands used in R for the Wishart distribution (in particular, check the help session for the command rWishart());
• some other commands and additional extensions are available in some packages which needs to be downloaded and installed separately, e.g. function Wishart() in the library ‘MCMCpack’, or Wishart() in the library ‘mixAK’, or rwishart() in the library ‘dlm’, dist.Wishart() in the library ‘LaplacesDemon’, and others;

A simple random sample from Wishart distribution $$W_{1}(\Sigma = 1, n = 10)$$ (equivalently a $$\chi^2$$ distribution with $$n = 10$$ degrees of freedom) can be obtained, for instance, as

S <- as.matrix(1)
sample <- rWishart(10, df = 10, S)

Using a standard approach for generating a univariate random sample we would use the R function rchisq() for this scenario. Try to use both and using some histograms compare whether do they correspond or not.For large sample size $$n \in \mathbb{N}$$ we they should be quite similar/identical.

set.seed(1234)
sampleWishart <- rWishart(5000, df = 10, S)
sampleChiSq <- rchisq(5000, df = 10)

par(mfrow = c(1,2))
hist(sampleWishart, col = "lightblue", main = expression(paste("Wishart Distribution", sep ="")), xlim = c(0, 40), breaks = 30, freq = F)
lines(density(sampleWishart), col = "red", lwd = 2)
hist(sampleChiSq, col = "lightblue", main = expression(paste(chi^2, "Distribution", sep = "")), xlim = c(0,40), breaks = 30, freq = F)
lines(density(sampleChiSq), col = "red", lwd = 2)

#### To Do by Yourself

• consider some Wishart distribution for some $$p \in \mathbb{N}$$ such that $$p > 1$$;
• generate a sample from the Wishart distribution $$W_{p}$$ with some reasonable parameters $$\Sigma$$ and $$df$$;
• consider some general matrix $$A \in \mathbb{R}^{p \times q}$$ for some $$q \in \mathbb{N}$$ such that $$p \neq q$$ and verify (e.g. by histograms) that $$A^\top \mathbb{M} A \sim W_{q}(A^\top\Sigma A, n)$$, where $$\mathbb{M} \sim W_{p}(\Sigma, n)$$;

#### Hotelling’s $$\boldsymbol{T^2}$$ Distribution

In a similar manner as we define a classical $$t$$-distribution in an univariate case (i.e. standard normal $$N(0,1)$$ variable devided by a square root of a $$\chi^2$$ variable normalized by its degrees of freedom) we define a multivariate generalization (Hotelling’s $$T^{2}$$ distribution) as $n \boldsymbol{Y}^{\top} \mathbb{M}^{-1} \boldsymbol{Y} \sim T^{2}(n, p),$ to be a random variable with the Hotelling’s $$T^2$$ distribution with $$p \in \mathbb{N}$$ to be the dimension of $$Y \sim N_{p}(0, \mathbb{I})$$ and $$n \in \mathbb{N}$$ being the parameter of the Wishart distribution of $$\mathbb{M} \sim W_{p}(\mathbb{I}, n)$$.

A special case for $$p = 1$$ gives the standard Fisher distribution with one and $$n$$ degrees of freedom (equivalently a square of the $$t$$-distribution with $$n$$ degrees of freedom). The Hotelling’s $$T^2$$ distribution with parameters $$p, n \in \mathbb{N}$$ can be identivied with the Fisher distribution using the following expression $T^{2}(p, n) \equiv \frac{n p}{n - p + 1}F_{p, n - p + 1}.$

In general, the relationship betwenn the Hotelling’s $$T^2$$ distribution and the Fisher distribution is given by the same expression $T^{2}(p, n) \equiv \frac{n p}{n - p + 1}F_{p, n - p + 1}.$

The effect of different parameter settings can be (a little) visualized, for instance, by the following figure:

samples <- NULL
samples <- cbind(samples, rf(n = 1000, df1 = 1, df2 = 1))
samples <- cbind(samples, rf(n = 1000, df1 = 1, df2 = 10))

samples <- cbind(samples, rf(n = 1000, df1 = 10, df2 = 10))
samples <- cbind(samples, rf(n = 1000, df1 = 10, df2 = 100))

samples <- cbind(samples, rf(n = 1000, df1 = 100, df2 = 10))
samples <- cbind(samples, rf(n = 1000, df1 = 100, df2 = 100))

samples <- cbind(samples, rf(n = 1000, df1 = 1000, df2 = 100))
samples <- cbind(samples, rf(n = 1000, df1 = 1000, df2 = 1000))

par(mfrow = c(4,2))
for (i in 1:dim(samples)[2]){
hist(samples[,i], xlab=paste("Sample no.", i, sep = ""), col = "yellow", main = "", freq = F, breaks = 20)
lines(density(samples[,i]), col = "red")
}

Moreover, using the relationship between the Hotelling’s $$T^2$$ distribution and the Fisher’s F distriubtion we can effectively use the $$F$$-distribution to draw the corresponding critical values when needed for some statistical tests.

For some multivariate random sample $$\boldsymbol{X}_{1}, \dots, \boldsymbol{X}_{n} \sim N_{p}(\boldsymbol{\mu}, \Sigma)$$ for some unknown mean vector $$\boldsymbol{\mu} \in \mathbb{R}^{p}$$ and some variance-covariance matrix $$\Sigma$$, we have that $(n - 1)\Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big) \sim T^2(p, n - 1),$ which can be also equivalently expressed as $\frac{n - p}{p} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big) \sim F_{p, n - p}.$ This can be now used to construct confidence regions for the unknown mean vector $$\boldsymbol{\mu}$$ of for testing hypothesis about the true value of the vector of parameters $$\boldsymbol{\mu} = (\mu_{1}, \dots, \mu_{p})^{\top}$$. However, rather than construction a confidence region for $$\boldsymbol{\mu}$$ (which can be impractical in higher dimensions for even slightly larger values of $$p$$) one focusses on construction confidence intervals for the elements of $$\boldsymbol{\mu}$$ such that the mutual coverage is under control (usually we require a simultaneous coverage of $$(1 - \alpha)\times 100~\%$$ for some small $$\alpha \in (0,1)$$).

For a hypothesis test

$H_{0}: \boldsymbol{\mu} = \boldsymbol{\mu}_{0} \in \mathbb{R}^{p}$

$H_{1}: \boldsymbol{\mu} \neq \boldsymbol{\mu}_{0} \in \mathbb{R}^{p}$

we can use the following test statistic:

$(n - 1)\Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big),$

which, under the null hypothesis, follows the $$T^2(p, n - 1)$$ distribution. Equivalently we also have that

$\frac{n - p}{p} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)$

follows the Fisher $$F$$-distribution with $$p$$ and $$n - p$$ degrees of freedom.

In the R software we can use the library install.packages("DescTools") - if it is installed in the R environment, we can load the library by the command

library(DescTools)

and we can use the function HotellingsT2Test() to perform the test on multivariate mean vector (use the help session to find out how the function works).

Using the same approach we can also construct the confidence elipsoid for $$\boldsymbol{\mu} \in \mathbb{R}^p$$ - it holds that $\frac{n - p}{p} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big)^\top \mathcal{S}^{-1} \Big(\overline{\boldsymbol{X}} - \boldsymbol{\mu}_{0}\Big) \sim F_{p, n - p},$ and therefore, the following set $\left\{\boldsymbol{\mu} \in \mathbb{R}^p;~ \Big( \overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big)^\top \mathcal{S}^{-1} \Big( \overline{\boldsymbol{X}} - \boldsymbol{\mu}\Big) \leq \frac{p}{n - p} F_{p, n- p}(1 - \alpha) \right\}$ is a confidence region at the confidence level of $$\alpha = 1 - \alpha$$ for the vector of parameters $$\boldsymbol{\mu} \in \mathbb{R}^p$$ - an interior of iso-distance ellipsoid in $$\mathbb{R}^p$$.

A brief example on how it works (example from the lecture notes):

library(mvtnorm)
s=matrix(c(1,-0.5,-0.5,1),2);x=seq(-3,3,by=0.015) ### variance-covariance matrix

contour(x,x,outer(x,x, function(x


text:
# Math Help - A few questions about integrals

1. ## A few questions about integrals

Hello!

I am having problem solving this one:

∫ -y/(x^2+y^2)dx
∫t* sqrt(1+4t^2) dt

Finding g(x) in the following equation:

(g'(x)-g(x))x sin xy = 0

2. Originally Posted by tinyone
Hello!

I am having problem solving this one:

∫ -y/(x^2+y^2)dx
Think of trigonometric functions. You'll have to use Inverse tangent function here..

Originally Posted by tinyone
∫t* sqrt(1+4t^2) dt
Use the substitution rule.. let $u=4t^2+1$

show your work if you are having any problems.


text:
## STOKE’S THEOREM

Let $S$ be an open, two-sided surface bounded by a closed non-intersecting curve $C$ (simple closed curve). Consider a directed line normal to $S$ as positive if it is on one side of $S$, and negative if it is on the other side of $S$. The choice of which side is positive is arbitrary but should be decided upon in advance. Call the direction or sense of $C$ positive if an observer, walking on the boundary of $S$ with his head pointing in the direction of the positive normal, has the surface on his left. Then if $A_1,\ A_2,\ A_3$ are single-valued, continuous, and have continuous first partial derivatives in a region of space including $S$, we have

$\displaystyle \int_C[A_1dx + A_2dy + A_3dz] =\\\vspace{0.2in} \underset{S}{\iint}\left[ \left( \frac{\partial A_3}{\partial y} -\frac{\partial A_2}{\partial z} \right)\cos\alpha + \left( \frac{\partial A_1}{\partial z} -\frac{\partial A_3}{\partial x} \right)\cos\beta + \left( \frac{\partial A_2}{\partial x} -\frac{\partial A_1}{\partial y} \right)\cos\gamma \right]dS \cdots(38)$

In vector form with $\bold{A} = A_1\bold{i} + A_2\bold{j} + A_3\bold{k}$ and $\bold{n} = \cos\alpha\bold{i} + \cos\beta\bold{j} + \cos\gamma\bold{k}$, this is simply expressed as

$\displaystyle \int_C \bold{A}\cdot d\bold{r} = \underset{S}{\iint}(\nabla\times\bold{A})\cdot\bold{n}dS\cdots(39)$

In words this theorem, called Stoke’s theorem, states that the line integral of the tangential component of a vector $\bold{A}$ taken around a simple closed curve $C$ is equal to the surface integral of the normal component of the curl of $A$ taken over any surface $S$ having $C$ as a boundary. Note that if, as a special case $\nabla\times\bold{A} = 0$ in (39), we obtain the result (28).

## EVALUATION OF LINE INTEGRALS

If the equation of a curve C in the plane $z = 0$ is given as $y = f(x)$, the line integral (14) is evaluated by placing $y = f(x),\ dy = f'(x)dx$ in the integrand to obtain the definite integral

$\displaystyle \int_{a_1}^{a_2}[P\{x, f(x)\}dx + Q\{x, f(x)\}f'(x)dx] \cdots(19)$

which is then evaluated in the usual manner.

Similarly if C is given as $x = g(y)$, then $dx = g'(y)dy$ and the line integral becomes

$\displaystyle \int_{b_1}^{b_2}[P\{g(y), y\}g'(y)dy + Q\{g(y), y\}dy]\cdots(20)$

If C is given in parametric form $x = \phi(t),\ y = \psi(t)$, the line integral becomes

$\displaystyle \int_{t_1}^{t_2} [P\{ \phi(t),\ \psi(t) \}\phi'(t)dt + Q\{ \phi(t),\ \psi(t) \}\psi'(t)dt] \cdots (21)$

where $t_1$ and $t_2$ denote the values of $t$ corresponding to points $A$ and $B$ respectively.

Combination of the above methods may be used in the evaluation.

Similar methods are used for evaluating line integrals along space curve.

## VECTOR NOTATION FOR LINE INTEGRALS

It is often convenient to express a line integral in vector form as an aid in physical or geometric understanding as well as for brevity of notation. For example, we can express the line integral (15) in the form

$\displaystyle \int_{C}[A_1dx + A_2dy + A_3dz] \\ = \int_{C} (A_1\bold{i} + A_2\bold{j} + A_3\bold{k}) \cdot (dx\bold{i} + dy\bold{j} + dz\bold{k})\\ = \int_{C} \bold{A}\cdot d\bold{r} \cdots (17)$

where $\bold{A} = A_1\bold{i} + A_2\bold{j} + A_3\bold{k}$ and $d\bold{r} = dx\bold{i} + dy\bold{j} + dz\bold{k}$. The line integral (14) is a special case of this with $z = 0$.

If at each point (x, y, z) we associate a force F acting on an object (i.e. if a force field is defined), then

$\displaystyle \int_{C} \bold{F}\cdot d\bold{r}\cdots(18)$

represents physically the total work done in moving the object along the curve C.

## LINE INTEGRALS

Let C be a curve in the xy plane which connects points $A (a_1, b_1)$ and $B (a_2, b_2)$, (see Fig. 6-2). Let $P(x, y)$ and $Q(x, y)$ be single-valued functions defined at all points of C. Subdivide C into n parts by choosing n – 1 points on it given by $(x_1, y_1),\ (x_2, y_2),\ \dots,\ (x_{n-1}, y_{n-1})$. Call $\Delta x_k = x_k - x_{k-1}$ and $\Delta y_k = y_k - y_{k-1},\ k = 1,\ 2,\ \dots\ n$ and suppose that points $(\xi_k, \eta_k)$ are chosen so that they are situated on C between points $(x_{k-1}, y_{k-1})$ and $(x_k, y_k)$. Form the sum

$\displaystyle \sum_{k=1}^{n}\{P(\xi_k, \eta_k)\Delta x_k + Q(\xi_k, \eta_k)\Delta y_k\}\cdots(13)$

The limit of this sum as $n\rightarrow\infty$ in such a way that all quantities $\Delta x_k,\ \Delta y$ approaches zero, if such limit exists, is called a line integral along C and is denoted by

$\displaystyle \int_C \left[ P(x, y)dx + Q(x, y)dy \right]$ or $\displaystyle \int_{(a_1, b_1)}^{(a_2, b_2)}\left[ Pdx + Qdy \right]\cdots(14)$

The limit does exist if P and Q are continuous (or piecewise continuous) at all points of C. The value of the integral depends in general on P, Q, the particular curve C, and on the limits $(a_1, b_1)$ and $(a_2, b_2)$.

In an exactly analogous manner one may define a line integral along a curve C in three dimensional space as

$\displaystyle \lim\limits_{n \rightarrow\infty}\sum_{k=1}^{n}\left\{ A_1(\xi_k, \eta_k, \zeta_k)\Delta x_k + A_2(\xi_k, \eta_k, \zeta_k)\Delta y_k + A_3(\xi_k, \eta_k, \zeta_k)\Delta z_k \right\} \\ = \int_C \left[ A_1dx + A_2dy + A_3dz \right] \cdots(15)$

where $A_1$, $A_2$ and $A_3$ are functions of $x$, $y$ and $z$.

Other types of line integrals, depending on particular curves, can be defined. For example, if $\Delta s_k$ denotes the arc length along curve C in the above figure between points $(x_k, y_k)$ and $(x_{k+1}, y_{k+1})$, then

$\displaystyle \lim\limits_{n \rightarrow \infty} \sum_{k=1}^{n} U(\xi_k, \eta_k)\Delta s_k = \int_C U(x, y)ds\cdots(16)$

is called the line integral of $U(x, y)$ along curve C. Extensions to three (or higher) dimensions are possible.


text:
Circle AMC10/12
2009

Problem - 1553
Circle $A$ has radius $100$. Circle $B$ has an integer radius $r<100$ and remains internally tangent to circle $A$ as it rolls once around the circumference of circle $A$. The two circles have the same points of tangency at the beginning and end of cirle $B$'s trip. How many possible values can $r$ have?


text:
# Particular Solution

Consider the nonhomogeneous differential equation ${{y}}^{{{\left({n}\right)}}}+{a}_{{{n}-{1}}}{\left({x}\right)}{{y}}^{{{\left({n}-{1}\right)}}}+\ldots+{a}_{{1}}{\left({x}\right)}{y}'+{a}_{{0}}{\left({x}\right)}{y}={g{{\left({x}\right)}}}$.

Recall from the section about linear independence that the solution of the corresponding homogeneous equation ${{y}}^{{{\left({n}\right)}}}+{a}_{{{n}-{1}}}{\left({x}\right)}{{y}}^{{{\left({n}-{1}\right)}}}+\ldots+{a}_{{1}}{\left({x}\right)}{y}'+{a}_{{0}}{\left({x}\right)}{y}={0}$ is given as ${y}_{{h}}={c}_{{1}}{y}_{{1}}{\left({x}\right)}+{c}_{{2}}{y}_{{2}}{\left({x}\right)}+\ldots+{c}_{{n}}{y}_{{n}}{\left({x}\right)}$.

If ${y}_{{p}}$ is the solution of the nonhomogeneous equation, the general solution of the nonhomogeneous equation is ${y}={y}_{{h}}+{y}_{{p}}$.

Proof.

Plug ${y}={y}_{{h}}+{y}_{{p}}$ into the differential equation:

${{\left({y}_{{h}}+{y}_{{p}}\right)}}^{{{\left({n}\right)}}}+{a}_{{{n}-{1}}}{\left({x}\right)}{{\left({y}_{{h}}+{y}_{{p}}\right)}}^{{{\left({n}-{1}\right)}}}+\ldots+{a}_{{1}}{\left({x}\right)}{\left({y}_{{h}}+{y}_{{p}}\right)}'+{a}_{{0}}{\left({x}\right)}{y}=$

$={\left({{y}_{{h}}^{{{\left({n}\right)}}}}+{a}_{{{n}-{1}}}{\left({x}\right)}{{y}_{{h}}^{{{\left({n}-{1}\right)}}}}+\ldots+{a}_{{1}}{\left({x}\right)}{y}_{{h}}'+{a}_{{0}}{\left({x}\right)}{y}_{{h}}\right)}+{\left({{y}_{{p}}^{{{\left({n}\right)}}}}+{a}_{{{n}-{1}}}{\left({x}\right)}{{y}_{{p}}^{{{\left({n}-{1}\right)}}}}+\ldots+{a}_{{1}}{\left({x}\right)}{y}_{{p}}'+{a}_{{0}}{\left({x}\right)}{y}_{{p}}\right)}=$

$={0}+{g{{\left({x}\right)}}}={g{{\left({x}\right)}}}$

So, indeed, ${y}={y}_{{p}}+{y}_{{h}}$ is the solution of the nonhomogeneous differential equation.

Now, let's consider an example.

Example. The solution of the homogeneous equation ${y}''+{y}'-{2}{y}={0}$ is ${y}_{{h}}={c}_{{1}}{{e}}^{{-{2}{x}}}+{c}_{{2}}{{e}}^{{{x}}}$. For the nonhomogeneous equation ${y}''+{y}'-{2}{y}={2}{\left({1}+{x}-{{x}}^{{2}}\right)}$, the particular solution is ${y}_{{p}}={{x}}^{{2}}$. So, the general solution of the nonhomogeneous differential equation is ${y}''+{y}'-{2}{y}={2}{\left({1}+{x}-{{x}}^{{2}}\right)}$ is ${y}={y}_{{h}}+{y}_{{p}}={c}_{{1}}{{e}}^{{-{2}{x}}}+{c}_{{2}}{{e}}^{{x}}+{{x}}^{{2}}$.


text:
If the angle 100^circ is in standard position, which quadrant does it terminate?

quadrant 2 because it is greater than ${90}^{\circ}$ and less than ${180}^{\circ}$