add `rand` for `ProbabilitySimplex` #604

adienes · 2023-05-10T23:00:46Z

uniformly samples from probability simplex if no point given. otherwise returns a mean-zero vector of gaussians.

codecov · 2023-05-10T23:21:12Z

Codecov Report

Merging #604 (5308add) into master (836e716) will increase coverage by 0.00%.
The diff coverage is 94.73%.

@@           Coverage Diff           @@
##           master     #604   +/-   ##
=======================================
  Coverage   99.01%   99.01%           
=======================================
  Files         104      106    +2     
  Lines       10138    10187   +49     
=======================================
+ Hits        10038    10087   +49     
  Misses        100      100

Impacted Files	Coverage Δ
src/Manifolds.jl	`86.66% <ø> (ø)`
src/manifolds/ProbabilitySimplexEuclideanMetric.jl	`90.00% <90.00%> (ø)`
src/manifolds/ProbabilitySimplex.jl	`100.00% <100.00%> (ø)`

... and 3 files with indirect coverage changes

📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more

adienes · 2023-05-11T01:43:17Z

sorry should have squashed all those first. had to turn off the random tangent test since isapprox has atol set to 0, and thus isapprox(_, 0) is testing for exact equality.

I think pr is ok now? let me know

kellertuer

Thanks for your PR – contributions are always happily welcome :)

I have two small remarks – and never mind several commits – for a PR that is always good to see also a bit the development, we squash-merge PRs anyways.

src/manifolds/ProbabilitySimplex.jl

test/manifolds/probability_simplex.jl

mateuszbaran · 2023-05-11T07:16:19Z

sorry should have squashed all those first. had to turn off the random tangent test since isapprox has atol set to 0, and thus isapprox(_, 0) is testing for exact equality.

The easiest way to handle this is setting rand_tvector_atol_multiplier=1 keyword argument to test_manifold, or a bigger number if needed. I think we should start setting non-zero absolute tolerance for tangent vector tests by default but that's a separate issue 🙂 .

kellertuer · 2023-05-11T07:17:51Z

I think we should start setting non-zero absolute tolerance for tangent vector tests by default but that's a separate issue 🙂 .

Good point to keep in mind for #548 .

src/manifolds/ProbabilitySimplex.jl

sethaxen · 2023-05-11T18:52:49Z

src/manifolds/ProbabilitySimplex.jl

+        _r = sort(vcat(rand(n - 1), [0.0, 1.0]))
+        pX .= (circshift(_r, -1) .- _r)[1:n]


We could make fewer allocations and avoid the O(n*log(n)) sorting by using the fact that Dirichlet(alpha=ones(n)) is the uniform distribution on the probability simplex, so the usual Gamma distribution (see https://en.wikipedia.org/wiki/Dirichlet_distribution#Random_variate_generation) approach for drawing from Dirichlet reduces to drawing from Exponential(1):

Suggested change

_r = sort(vcat(rand(n - 1), [0.0, 1.0]))

pX .= (circshift(_r, -1) .- _r)[1:n]

Random.randexp!(pX)

LinearAlgebra.normalize!(pX, 1)

Nice! That should also be reflected in the docs.

sethaxen · 2023-05-11T18:55:36Z

src/manifolds/ProbabilitySimplex.jl

+        _r = randn(n) .* σ
+        pX .= _r .- mean(_r)


We don't need the allocation to r_.

Suggested change

_r = randn(n) .* σ

pX .= _r .- mean(_r)

Random.randn!(pX)

pX .= (pX .- mean(pX)) .* σ

sethaxen · 2023-05-11T18:55:48Z

src/manifolds/ProbabilitySimplex.jl

+        _r = randn(n) .* σ
+        pX .= _r .- mean(_r)


We don't need the allocation to r_.

Suggested change

_r = randn(n) .* σ

pX .= _r .- mean(_r)

Random.randn!(pX)

pX .= (pX .- mean(pX)) .* σ

sethaxen · 2023-05-11T18:57:23Z

src/manifolds/ProbabilitySimplex.jl

+    end
+    return pX
+end
+function Random.rand!(


I know elsewhere in this package we write separate methods for with and without rng, but this was a mistake. In general, we should implement sampling methods with rng and then the one without rng should just dispatch to this by passing Random.default_rng() as the first argument. This is what Random itself does.

I was not aware f that we should then slowly improve our code :)

I will make a PR that improves this.

kellertuer · 2023-05-12T08:07:11Z

With your current tolerance factor of 1 you are still checking for machine precision zero, which might be a bit tough (that is why 1.6 failed, maybe raise the factor to 10.0?

test/manifolds/probability_simplex.jl

kellertuer · 2023-05-12T12:50:32Z

From my side this nearly looks good now, just that the doc string could maybe include a math formula or a reference where to find more information (maybe Seth's wikipedia remark)?

sethaxen · 2023-05-12T13:18:39Z

There's a reference here I haven't read: https://stackoverflow.com/a/67202070

sethaxen · 2023-05-12T15:44:47Z

Okay, got it. The standard reference is Theorem 2.2, Chapter 5 of

Luc Devroye. Non-Uniform Random Variate Generation. Springer-Verlag, New York, 1986. doi: 10.1007/978-1-4613-8643-8

kellertuer · 2023-05-12T16:36:53Z

Cool, maybe adding that as a footnote reference would be cool (at some point we should aim for some bib package in the docs, let's hope some package like that soon magically appears :D )

src/manifolds/ProbabilitySimplex.jl

adienes · 2023-05-12T17:02:38Z

thank you for adding the reference. we should probably merge one of this PR or #605 first, then rebase the other? I think as-is they are slightly incompatible

kellertuer · 2023-05-12T17:19:24Z

Oh I might have crashed formatting, sorry for that (editing here has its disadvantages). But I‘ll activate Docs rendering as well now :)

There seems to be a space too much in the footnote (i.e. [^ D...] (the one after ^); and if you could run formatter locally that would be great: in the Manifolds.jl folder run using JuliaFormatter; format(".")

mateuszbaran · 2023-05-12T17:30:29Z

thank you for adding the reference. we should probably merge one of this PR or #605 first, then rebase the other? I think as-is they are slightly incompatible

I think we can merge this one first, and then I will fix incompatibility in #605 .

src/manifolds/ProbabilitySimplex.jl

kellertuer

My points are clarified. Thanks for the contribution and the good further ideas.

sethaxen

Just a minor docstring fix. Also, this should get a patch version bump. Otherwise, LGTM.

src/manifolds/ProbabilitySimplex.jl

mateuszbaran · 2023-05-16T17:26:56Z

I think I will merge this and bump version in #605 -- thank you for your contribution 🙂

kellertuer · 2023-05-16T17:29:08Z

But we already discussed, that this is w.r.t. the „wrong” metric, i.e. not our default one but the Euclidean (with respect to which the manifold is not complete).

mateuszbaran · 2023-05-16T17:40:15Z

Yes, it's not uniform wrt. our metric but we never claimed to always use uniform distributions (for manifolds of finite volume). And we can change it later.

kellertuer · 2023-05-16T18:18:53Z

Sure, I would still prefer if we could switch this one to be the one with the metric manifold / euclidean metric?

mateuszbaran · 2023-06-13T13:35:32Z

Wait... that's actually a good question. I assumed change_representer and change_metric are fine but when I look at them now they indeed seem to produce vectors that are not tangent. @kellertuer do you remember where you got the formulas from?

adienes · 2023-06-13T13:49:20Z

that commit is just a rebase on master, diff shouldn't have changed

kellertuer · 2023-06-13T13:51:51Z

We usually merge, that is easier (for PRs we squash-merge, which makes the master branch relatively clean again).

I will check the change of metric soon.

kellertuer · 2023-06-13T14:23:27Z

Hm, I am now a bit confused, since the documentation of exp for example follows https://arxiv.org/pdf/1603.05285.pdf, Proposition 2.4, but the implementation – at least as far as I understand it – does something completely different.

The metric conversion seems wrong for now, I might have stayed too close to what is done on the positive numbers (element wise), will check for a correction. The change of representer seems right, though.

edit: Ah no, the change representer is also wrong, but Proposition (2.3) provides the right formula, but then I can probably derive the change metric from that as well. - eh, oh! That formula includes projection which we do before, added that note to the docs, but change representer is correct. Will check how to adapt that though.

kellertuer · 2023-06-13T15:15:50Z

HM, that's a tough one to solve, because it seems there are only n degrees of freedom (for the change metric) but I currently seem to have n+1 conditions.
Yes, I am convinced we have to remove that method because it is just not possible.

I have to dismiss my approval since there still seems to be work to be done

mateuszbaran · 2023-06-13T16:47:47Z

I've discussed with @kellertuer what is wrong with change_representer and change_metric. I will make a PR fixing them and then you will be able to continue 🙂 . Thank you for your question about validity of change_metric, I wouldn't notice that something is wrong without it.

kellertuer · 2023-06-13T18:36:44Z

Yeah, thanks for the discussion it seems there might be a change metric :)

mateuszbaran · 2023-06-14T09:45:50Z

change_metric is now fixed on master branch, so you can continue now 🙂

…lex (JuliaManifolds#627)

* Add a few exports, mainly reexports from ManifoldsBase. * Improve caching. i.e. doing it correctly. * 📖🩹 at least a bit, * 📖🩹 2 * export FisherRaoMetric.

adienes · 2023-06-17T01:02:48Z

I had tried also to add some tests for the ProbabilitySimplexEuclideanMetric, but it was very hard to get any to pass since basically only rand is defined. it seems to want retract ...

mateuszbaran · 2023-06-17T07:18:46Z

Hm, I think you can just take exp from Euclidean:

exp!(::MetricManifold{ℝ,<:ProbabilitySimplex,<:EuclideanMetric}, q, p, X) = (q .= p .+ X)
exp!(::MetricManifold{ℝ,<:ProbabilitySimplex,<:EuclideanMetric}, q, p, X, t::Number) = (q .= p .+ t .* X)

and it should more or less work.

adienes · 2023-06-17T15:11:13Z

g2g?

test/manifolds/probability_simplex.jl

src/manifolds/ProbabilitySimplexEuclideanMetric.jl

src/manifolds/ProbabilitySimplex.jl

Co-authored-by: Seth Axen <seth@sethaxen.com>

kellertuer · 2023-06-24T22:28:40Z

But then we should also register 0.8.71 soon, because I think we were just waiting for this PR.

mateuszbaran · 2023-06-25T12:48:13Z

Yes, I've even forgotten we have a few changes on master without a new registered version.

adienes force-pushed the randpr branch from 8b5c167 to 8074a7b Compare May 11, 2023 01:41

kellertuer requested changes May 11, 2023

View reviewed changes

src/manifolds/ProbabilitySimplex.jl Outdated Show resolved Hide resolved

test/manifolds/probability_simplex.jl Outdated Show resolved Hide resolved

kellertuer reviewed May 11, 2023

View reviewed changes

src/manifolds/ProbabilitySimplex.jl Outdated Show resolved Hide resolved

sethaxen requested changes May 11, 2023

View reviewed changes

mateuszbaran mentioned this pull request May 11, 2023

Simplify random point/tangent generation by using default_rng #605

Merged

kellertuer reviewed May 12, 2023

View reviewed changes

test/manifolds/probability_simplex.jl Outdated Show resolved Hide resolved

test/manifolds/probability_simplex.jl Outdated Show resolved Hide resolved

kellertuer reviewed May 12, 2023

View reviewed changes

src/manifolds/ProbabilitySimplex.jl Show resolved Hide resolved

kellertuer added the preview docs Add this label if you want to see a PR-preview of the documentation label May 12, 2023

kellertuer reviewed May 12, 2023

View reviewed changes

src/manifolds/ProbabilitySimplex.jl Outdated Show resolved Hide resolved

kellertuer previously approved these changes May 12, 2023

View reviewed changes

mateuszbaran approved these changes May 12, 2023

View reviewed changes

kellertuer requested a review from sethaxen May 13, 2023 19:46

sethaxen requested changes May 13, 2023

View reviewed changes

src/manifolds/ProbabilitySimplex.jl Outdated Show resolved Hide resolved

src/manifolds/ProbabilitySimplex.jl Outdated Show resolved Hide resolved

sethaxen mentioned this pull request May 14, 2023

parallel_transport_to for ProbabilitySimplex #606

Closed

type params simpler

26592dc

adienes force-pushed the randpr branch from 248d166 to 26592dc Compare June 13, 2023 13:48

mateuszbaran mentioned this pull request Jun 13, 2023

Fix change_metric, change_representer and project for ProbabilitySimplex #627

Merged

mateuszbaran and others added 4 commits June 16, 2023 20:52

Fix change_metric, change_representer and project for ProbabilitySimp…

0e314c6

…lex (JuliaManifolds#627)

Add Exports, Improve GHA Caches. (JuliaManifolds#626)

9f00db7

* Add a few exports, mainly reexports from ManifoldsBase. * Improve caching. i.e. doing it correctly. * 📖🩹 at least a bit, * 📖🩹 2 * export FisherRaoMetric.

new change metric

898c841

Merge branch 'master' into randpr

d6594ee

add exp and test

6dfe32a

mateuszbaran reviewed Jun 17, 2023

View reviewed changes

test/manifolds/probability_simplex.jl Outdated Show resolved Hide resolved

Update test/manifolds/probability_simplex.jl

9bc3bca

mateuszbaran reviewed Jun 17, 2023

View reviewed changes

src/manifolds/ProbabilitySimplexEuclideanMetric.jl Outdated Show resolved Hide resolved

Update src/manifolds/ProbabilitySimplexEuclideanMetric.jl

490d3de

mateuszbaran approved these changes Jun 17, 2023

View reviewed changes

sethaxen reviewed Jun 19, 2023

View reviewed changes

src/manifolds/ProbabilitySimplex.jl Outdated Show resolved Hide resolved

Update src/manifolds/ProbabilitySimplex.jl

5308add

Co-authored-by: Seth Axen <seth@sethaxen.com>

mateuszbaran merged commit d24f6d2 into JuliaManifolds:master Jun 24, 2023
25 of 30 checks passed

		_r = sort(vcat(rand(n - 1), [0.0, 1.0]))
		pX .= (circshift(_r, -1) .- _r)[1:n]

add rand for ProbabilitySimplex #604

add rand for ProbabilitySimplex #604

Conversation

adienes commented May 10, 2023

codecov bot commented May 10, 2023 • edited

Codecov Report

adienes commented May 11, 2023

kellertuer left a comment

Choose a reason for hiding this comment

mateuszbaran commented May 11, 2023

kellertuer commented May 11, 2023

sethaxen May 11, 2023

Choose a reason for hiding this comment

kellertuer May 11, 2023

Choose a reason for hiding this comment

sethaxen May 11, 2023

Choose a reason for hiding this comment

sethaxen May 11, 2023

Choose a reason for hiding this comment

sethaxen May 11, 2023

Choose a reason for hiding this comment

kellertuer May 11, 2023

Choose a reason for hiding this comment

mateuszbaran May 11, 2023

Choose a reason for hiding this comment

kellertuer commented May 12, 2023

kellertuer commented May 12, 2023

sethaxen commented May 12, 2023

sethaxen commented May 12, 2023 • edited

kellertuer commented May 12, 2023

adienes commented May 12, 2023

kellertuer commented May 12, 2023 • edited

mateuszbaran commented May 12, 2023

kellertuer left a comment

Choose a reason for hiding this comment

sethaxen left a comment

Choose a reason for hiding this comment

mateuszbaran commented May 16, 2023

kellertuer commented May 16, 2023

mateuszbaran commented May 16, 2023

kellertuer commented May 16, 2023

mateuszbaran commented Jun 13, 2023

adienes commented Jun 13, 2023

kellertuer commented Jun 13, 2023

kellertuer commented Jun 13, 2023 • edited

kellertuer commented Jun 13, 2023 • edited

mateuszbaran commented Jun 13, 2023

kellertuer commented Jun 13, 2023

mateuszbaran commented Jun 14, 2023

adienes commented Jun 17, 2023

mateuszbaran commented Jun 17, 2023

adienes commented Jun 17, 2023

kellertuer commented Jun 24, 2023

mateuszbaran commented Jun 25, 2023

add `rand` for `ProbabilitySimplex` #604

add `rand` for `ProbabilitySimplex` #604

codecov bot commented May 10, 2023 •

edited

sethaxen commented May 12, 2023 •

edited

kellertuer commented May 12, 2023 •

edited

kellertuer commented Jun 13, 2023 •

edited

kellertuer commented Jun 13, 2023 •

edited