Skip to content

Commit

Permalink
I work out where the sqrt(2) comes from in the Langevin equation.
Browse files Browse the repository at this point in the history
  • Loading branch information
pfh committed Mar 13, 2022
1 parent ad73209 commit 38e7ecb
Show file tree
Hide file tree
Showing 3 changed files with 17 additions and 8 deletions.
10 changes: 7 additions & 3 deletions inst/htmlwidgets/lib/langevitour.js
Original file line number Diff line number Diff line change
Expand Up @@ -957,7 +957,7 @@ class Langevitour {

compute(realElapsed) {
let damping = 0.2 *Math.pow(10, this.get('dampInput').value);
let temperature = 0.05 *Math.pow(10, this.get('tempInput').value);
let temperature = 0.1 *Math.pow(10, this.get('tempInput').value);
let repulsion = 1.0 *Math.pow(10, this.get('repulsionInput').value);
let attraction = 1.0 *Math.pow(10, this.get('labelInput').value);
let doTemp = this.get('tempCheckbox').checked;
Expand Down Expand Up @@ -987,10 +987,14 @@ class Langevitour {

if (doTemp) {
// Damping reduces the variance * velKeep^2
// We need to add velReplaceVar * desired steady state variance of temperature*2
// We need to add velReplaceVar * desired steady state variance of temperature

// Note the sqrt(2) that appears in the continuous form of the Langevin equation arises from
// sqrt(1-velKeep*velKeep) approx= sqrt(2*elapsed*damping) for small elapsed*damping

let velReplaceVar = 1 - velKeep*velKeep;
let noise = times(proj.length, times, this.m,
jStat.normal.sample, 0, Math.sqrt(2*temperature*velReplaceVar));
jStat.normal.sample, 0, Math.sqrt(temperature*velReplaceVar));

noise = removeSpin(noise, proj);

Expand Down
10 changes: 7 additions & 3 deletions langevitour.js
Original file line number Diff line number Diff line change
Expand Up @@ -957,7 +957,7 @@ class Langevitour {

compute(realElapsed) {
let damping = 0.2 *Math.pow(10, this.get('dampInput').value);
let temperature = 0.05 *Math.pow(10, this.get('tempInput').value);
let temperature = 0.1 *Math.pow(10, this.get('tempInput').value);
let repulsion = 1.0 *Math.pow(10, this.get('repulsionInput').value);
let attraction = 1.0 *Math.pow(10, this.get('labelInput').value);
let doTemp = this.get('tempCheckbox').checked;
Expand Down Expand Up @@ -987,10 +987,14 @@ class Langevitour {

if (doTemp) {
// Damping reduces the variance * velKeep^2
// We need to add velReplaceVar * desired steady state variance of temperature*2
// We need to add velReplaceVar * desired steady state variance of temperature

// Note the sqrt(2) that appears in the continuous form of the Langevin equation arises from
// sqrt(1-velKeep*velKeep) approx= sqrt(2*elapsed*damping) for small elapsed*damping

let velReplaceVar = 1 - velKeep*velKeep;
let noise = times(proj.length, times, this.m,
jStat.normal.sample, 0, Math.sqrt(2*temperature*velReplaceVar));
jStat.normal.sample, 0, Math.sqrt(temperature*velReplaceVar));

noise = removeSpin(noise, proj);

Expand Down
5 changes: 3 additions & 2 deletions vignettes/examples.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -86,7 +86,7 @@ A dataset from the [tourr](http://ggobi.github.io/tourr/) package.

By default, langevitour does not scale variables individually. This dataset rather requires scaling, so we supply our desired scaling.

I also demonstrate adding extra axes representing principal components and rotated principal components. The package GPArotation provides many factor rotation methods. The Bentler method seems to work well.
I also demonstrate adding extra axes representing principal components and rotated principal components. The package GPArotation provides many rotation methods. Using the Bentler method to produce sparse scores seems to highlight interesting structure.

```{r}
olives <- tourr::olive[,3:10]
Expand All @@ -95,13 +95,14 @@ olivesGroup <- paste(tourr::olive$region, tourr::olive$area)
# Find principal components
pca <- prcomp(sweep(olives, 2, olivesScale, "/"), rank=4)
olivesPC = sweep(pca$rotation, 1, olivesScale, "/")
olivesPC <- sweep(pca$rotation, 1, olivesScale, "/")
# Find a good rotation of the loadings from PCA
rotation <- GPArotation::bentlerT(pca$x)$Th
olivesR <- olivesPC %*% rotation
colnames(olivesR) <- paste0("R",seq_len(ncol(olivesR)))
langevitour(
olives,
group=olivesGroup,
Expand Down

0 comments on commit 38e7ecb

Please sign in to comment.