text:
### A Pluto.jl notebook ###
# v0.16.1

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : missing
        el
    end
end

# ╔═╡ b6681ace-5cb1-49cb-9ed7-e749c7a67428
using Distributions, Plots, StatsPlots, PlutoUI

# ╔═╡ a4a73fb0-2213-11ec-02f8-1d2409575d4c
md""" # Bayesian modelling

Three stages to building a model:
1. Data story: Motivate the model by narrating how the data might arise.
2. Update: Educate your model by feeding it the data.
3. Evaluate: evaluate the model, possibly leading to model revision.

### Example: Estimate the proportion of water on earth

Suppose you have a toy model of the Earth. You estimate the proportion of water on the model of Earth (and, if the model is a good one, the proportion of water on the actual Earth) by tossing it into the air and catching it. When you catch it record whether or not your right index finger is on land or water. Repeat this procedure.

###### 1. Data story
This story may either be *descriptive* specifying associations that can be used to predict outcomes, or *causal*, a theory of how some events produce other events.

For Earth example:
1. The true proportion of the Earth coverd in water is $P$.
2. A single toss and catch has a probability $p$ of producing a water (W) observation
3. Each toss of the globe is independent of each other.

###### 2. Update
Use the data to update the model.

For Earth example:
1. Assume all values of $P$ equally likely (the *prior*)
2. As observed data is binary, use Binomial distribution to find probability of observed data for each possible value of $P$
3. Update prior probability using Bayes Theorem and the observed data (the *posterior*)

Calculating these probabilities from the PDF's to generate the posterior is almost impossible for more complex examples, so numerical techniques used instead.

###### 3. Evaluate
The Bayesian model learns in a way that is demonstrably optimal. However this is only the performance of the model (the "small world") if the assumptions of the model do not match reality then the "large world" performance will not be good.

The objective is to check the model's adequacy for you purpose e.g. prediction, comprehension, measurement.

"""

# ╔═╡ 4dcf2480-7a37-4236-b51b-a171407bd8a0
md"""
#### Constructing the model
"""

# ╔═╡ e044742f-15dd-42e1-9271-43f11c7e2785
begin
	N = 9 # number of trials
	S = 6 # number of success (water)
end

# ╔═╡ c2d41bee-4952-429b-b3d8-aea46eefd2c4
md"""
Given the outcome of our experiment can be one of two options (W or L) with a fixed probability $p$ and each toss is independent (these come from the assumptions/data story), the probability of a single toss being W can be described by a binomial distribution:

``W\sim Binomial(N,p)``

"""

# ╔═╡ ef92132a-211d-4e8f-8ab2-12ae2c4a9c49
md"""
Initially we don't know what the parameter $P$ is, so assume all values equally likely, i.e. a uniform distribution between 0 and 1 (0% to 100% covered in water):

``P\sim Uniform(0,1)``
"""

# ╔═╡ e2b04571-e4e5-49fe-823d-af0404eceab7
P = Uniform(0,1);

# ╔═╡ e2bb3b47-1973-4d40-91f2-a703470d1a6a
begin
	plot(P, lw=2, legend=:none, title="prior distribution of P", ylims=(0,1.2))
	#plot!([0.0 1.0; 0.0 1.0],[0.0 0.0; 1.0 1.0], color=1)
end

# ╔═╡ 8f3de7b2-2fa3-4aa3-ba81-243aa4722a3a
md"""#### Update the model with data

The model gets updated (the posterior is generated) using Bayes theorem:

``Pr(p\mid W,L) = \frac{Pr(W,L\mid p)Pr(p)}{Pr(W,L)}``

Analytical this would give a posterior probability density function, but using grid approximation the values given are probabilities (masses):
"""

# ╔═╡ d142a6c3-e11c-44b5-930f-16deae98e74f
@bind n Slider(5:100, default=100) # fine-ness of the grid

# ╔═╡ 1a193ab0-10b0-447d-a68a-0fc978c6e804
md"""n = $n"""

# ╔═╡ 33fc930d-8693-41f9-8614-dc71870bd56d
begin
	grid = range(0,1,length=n) # possible values for P
	prior = @. pdf(Uniform(0,1), grid) # prior
	likelihood = @. pdf(Binomial(N,grid), S) # probability of observed data for diff P
	posterior = likelihood .* prior # unnormalised posterior distribution
	posterior /= sum(posterior) 	# normalise posterior to sum to one -> probability
end;

# ╔═╡ 4d10cc40-213e-4c4e-bea8-ea318aee18f7
begin
	plot(grid, posterior, markershape=:circle, mc=:white, msc=1, legend=:none, 		title="posterior distribution", annotations=(0.8,0.5," $(round(grid[argmax(posterior)], digits=3))"), ylims=(0,0.6))
	vline!([grid[argmax(posterior)]], color=1, ls=:dash)
end

# ╔═╡ 9149e9e5-b0b9-4c22-8795-d982a3cacf6a
md"""##### Collecting new data"""

# ╔═╡ 361b07cc-6554-40d1-9448-6f4e4cd4bb71
@bind Nnew Slider(1:100, default=9)

# ╔═╡ af901869-37bf-4992-bfc3-ee8477a37b34
md"""Number of tosses: $Nnew"""

# ╔═╡ e82146d1-ea29-40a3-b89f-5d178f4a900a
@bind Snew Slider(0:Nnew, default=6)

# ╔═╡ 2161be7c-571e-46c9-b3c2-0543d5e3296f
md"""Number of tosses on water $Snew"""

# ╔═╡ da10276c-1381-4edb-9171-1ec0820a1809
let
	n = 500
	grid = range(0,1,length=n) # possible values for P
	prior = repeat([1],n) # relative priors
	likelihood = @. pdf(Binomial(Nnew,grid), Snew) # probability of S successes
	posterior = likelihood .* prior |> x -> x/sum(x) # posterior distribution
	plot(grid, posterior, lw=2, legend=:none, title="posterior distribution", ylims=(0,0.04))
end

# ╔═╡ b8e49bec-d582-4695-b7ac-6b0ccbca5e08
md"""##### Trying out different priors (using original data)"""

# ╔═╡ 2cccdace-62df-4771-aaa4-1dd4342335e4
md"""###### Step prior
```prior = [zeros(n÷2)..., ones(n÷2)...]```
"""

# ╔═╡ a673dd91-ccdf-4783-8b24-357a7307d7d3
let
	n = 500
	grid = range(0,1,length=n) # possible values for P
	prior = [zeros(n÷2)..., ones(n÷2)...]
	likelihood = @. pdf(Binomial(Nnew,grid), Snew) # probability of S successes
	posterior = likelihood .* prior |> x -> x/sum(x) # posterior distribution
	p1 = plot(grid, prior, lw=2, legend=:none, title="prior distribution")
	p2 = plot(grid, posterior, lw=2, legend=:none, title="posterior distribution", ylims=(0,0.04))
	plot(p1,p2, layout=(1,2))
end

# ╔═╡ 3a34068e-4fed-40c3-89ef-2aabafb5144e
md"""###### Exponential prior
```prior = exp.(-5*abs.(collect(grid).-0.5) ) ```
"""

# ╔═╡ 9c599a59-96f7-4b2c-b690-13750a14eef8
let
	n = 500
	grid = range(0,1,length=n) # possible values for P
	prior = exp.(-5*abs.(collect(grid).-0.5) )
	likelihood = @. pdf(Binomial(Nnew,grid), Snew) # probability of S successes
	posterior = likelihood .* prior |> x -> x/sum(x) # posterior distribution
	p1 = plot(grid, prior, lw=2, legend=:none, title="prior distribution")
	p2 = plot(grid, posterior, lw=2, legend=:none, title="posterior distribution", ylims=(0,0.04))
	plot(p1,p2, layout=(1,2))
end

# ╔═╡ ffc34a30-9dba-40d0-afca-2279f4753661
md"""### Practice problems"""

# ╔═╡ a70c520f-860f-48be-a1a6-e87607298c37
md""" #### 2H1

##### Question
Suppose there are two species of panda, identical in every way except family sizes. Species A gives birth to twins 10% of the time, otherwise a single infant. Species B gives birth to twins 20% of the time, otherwise a single infant.

You have a new panda of unknown species just given birth to twins, what is the probability the next birth will also be twins?

##### Solution

1. Construct data story

Need to find probability of new panda being species A or B. The probability of next birth being twins will then be:

``Pr(A)*Pr(twins \mid A) + Pr(B)*Pr(twins \mid B)``

Species distribution:
- Assume prior to be equally likely for both species.
- Model outcome with binomial distribution: ``T \sim Binomial(2,p)`` where ``p = 0.1`` if species A or ``p = 0.2`` if species B.
"""

# ╔═╡ 8c493e53-bb7b-4110-a732-1cf93ecb754b
let
	N = 1
	S = 1
	grid = [0.1, 0.2]
	prior = [1, 1]
	likelihood = @. pdf(Binomial(N,grid), S)
	posterior = posterior = likelihood .* prior |> x -> x/sum(x)
	p = posterior .* grid |> sum
	md"""**Probability of next birth being twins = $p**"""
end

# ╔═╡ 9908c97f-efb4-4634-aa97-9b4d2f9d7acc
md""" #### 2H2

##### Question

Same as above, compute probability that the new panda is of species A given that the first birth was twins.

##### Answer

Want to find ``Pr(A \mid twins)`` which is reverse of ``Pr(twins \mid A)`` which has probability 0.1. Therefore use Bayes Theorem:

``Pr(A \mid twins) = \frac{Pr(twins \mid A)*Pr(A)}{Pr(twins)}``

"""

# ╔═╡ d690e1d2-654c-4da4-9255-d3871429d262
let
	p_twins_A = 0.1
	p_A = 0.5 # assuming equally likely
	p_twins = p_A*0.1 + (1-p_A)*0.2
	p_A_twins = (p_twins_A * p_A) / p_twins
	md"""**Probability of being species A given twins = $(p_A_twins)**"""
end

# ╔═╡ ee5bfc95-bbb0-4d1b-8a83-77e884171f05
md""" #### 2H3

##### Question

Same as above, the new panda next gives birth to a single infant. Compute probability that she is of species A.

##### Solution

Use posterior distribution from 2H1 which gives probabilities after first birth of twins as the prior and generate new posterior given new birth of not twins.

"""

# ╔═╡ b7033638-277e-4610-89ef-3abe33ef35bf
let
	N = 1
	S = 0
	grid = [0.1, 0.2]
	prior = [1, 2] # from 2H1
	likelihood = @. pdf(Binomial(N, grid), S) # probability of observing data given p
	posterior = likelihood .* prior |> x -> x ./ sum(x)
	md"""**Probability of being species A given second birth of twins = $(posterior[1])**"""
end

# ╔═╡ 00000000-0000-0000-0000-000000000001
PLUTO_PROJECT_TOML_CONTENTS = """
[deps]
Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
PlutoUI = "7f904dfe-b85e-4ff6-b463-dae2292396a8"
StatsPlots = "f3b207a7-027a-5e70-b257-86293d7955fd"

[compat]
Distributions = "~0.25.17"
Plots = "~1.22.3"
PlutoUI = "~0.7.14"
StatsPlots = "~0.14.28"
"""

# ╔═╡ 00000000-0000-0000-0000-000000000002
PLUTO_MANIFEST_TOML_CONTENTS = """
# This file is machine-generated - editing it directly is not advised

[[AbstractFFTs]]
deps = ["LinearAlgebra"]
git-tree-sha1 = "485ee0867925449198280d4af84bdb46a2a404d0"
uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c"
version = "1.0.1"

[[Adapt]]
deps = ["LinearAlgebra"]
git-tree-sha1 = "84918055d15b3114ede17ac6a7182f68870c16f7"
uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
version = "3.3.1"

[[ArgTools]]
uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f"

[[Arpack]]
deps = ["Arpack_jll", "Libdl", "LinearAlgebra"]
git-tree-sha1 = "2ff92b71ba1747c5fdd541f8fc87736d82f40ec9"
uuid = "7d9fca2a-8960-54d3-9f78-7d1dccf2cb97"
version = "0.4.0"

[[Arpack_jll]]
deps = ["Libdl", "OpenBLAS_jll", "Pkg"]
git-tree-sha1 = "e214a9b9bd1b4e1b4f15b22c0994862b66af7ff7"
uuid = "68821587-b530-5797-8361-c406ea357684"
version = "3.5.0+3"

[[Artifacts]]
uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33"

[[AxisAlgorithms]]
deps = ["LinearAlgebra", "Random", "SparseArrays", "WoodburyMatrices"]
git-tree-sha1 = "a4d07a1c313392a77042855df46c5f534076fab9"
uuid = "13072b0f-2c55-5437-9ae7-d433b7a33950"
version = "1.0.0"

[[Base64]]
uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"

[[Bzip2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "19a35467a82e236ff51bc17a3a44b69ef35185a2"
uuid = "6e34b625-4abd-537c-b88f-471c36dfa7a0"
version = "1.0.8+0"

[[Cairo_jll]]
deps = ["Artifacts", "Bzip2_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "JLLWrappers", "LZO_jll", "Libdl", "Pixman_jll", "Pkg", "Xorg_libXext_jll", "Xorg_libXrender_jll", "Zlib_jll", "libpng_jll"]
git-tree-sha1 = "f2202b55d816427cd385a9a4f3ffb226bee80f99"
uuid = "83423d85-b0ee-5818-9007-b63ccbeb887a"
version = "1.16.1+0"

[[ChainRulesCore]]
deps = ["Compat", "LinearAlgebra", "SparseArrays"]
git-tree-sha1 = "1417269aa4238b85967827f11f3e0ce5722b7bf0"
uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"
version = "1.7.1"

[[Clustering]]
deps = ["Distances", "LinearAlgebra", "NearestNeighbors", "Printf", "SparseArrays", "Statistics", "StatsBase"]
git-tree-sha1 = "75479b7df4167267d75294d14b58244695beb2ac"
uuid = "aaaa29a8-35af-508c-8bc3-b662a17a0fe5"
version = "0.14.2"

[[ColorSchemes]]
deps = ["ColorTypes", "Colors", "FixedPointNumbers", "Random"]
git-tree-sha1 = "a851fec56cb73cfdf43762999ec72eff5b86882a"
uuid = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
version = "3.15.0"

[[ColorTypes]]
deps = ["FixedPointNumbers", "Random"]
git-tree-sha1 = "024fe24d83e4a5bf5fc80501a314ce0d1aa35597"
uuid = "3da002f7-5984-5a60-b8a6-cbb66c0b333f"
version = "0.11.0"

[[Colors]]
deps = ["ColorTypes", "FixedPointNumbers", "Reexport"]
git-tree-sha1 = "417b0ed7b8b838aa6ca0a87aadf1bb9eb111ce40"
uuid = "5ae59095-9a9b-59fe-a467-6f913c188581"
version = "0.12.8"

[[Compat]]
deps = ["Base64", "Dates", "DelimitedFiles", "Distributed", "InteractiveUtils", "LibGit2", "Libdl", "LinearAlgebra", "Markdown", "Mmap", "Pkg", "Printf", "REPL", "Random", "SHA", "Serialization", "SharedArrays", "Sockets", "SparseArrays", "Statistics", "Test", "UUIDs", "Unicode"]
git-tree-sha1 = "31d0151f5716b655421d9d75b7fa74cc4e744df2"
uuid = "34da2185-b29b-5c13-b0c7-acf172513d20"
version = "3.39.0"

[[CompilerSupportLibraries_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae"

[[Contour]]
deps = ["StaticArrays"]
git-tree-sha1 = "9f02045d934dc030edad45944ea80dbd1f0ebea7"
uuid = "d38c429a-6771-53c6-b99e-75d170b6e991"
version = "0.5.7"

[[DataAPI]]
git-tree-sha1 = "cc70b17275652eb47bc9e5f81635981f13cea5c8"
uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a"
version = "1.9.0"

[[DataStructures]]
deps = ["Compat", "InteractiveUtils", "OrderedCollections"]
git-tree-sha1 = "7d9d316f04214f7efdbb6398d545446e246eff02"
uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
version = "0.18.10"

[[DataValueInterfaces]]
git-tree-sha1 = "bfc1187b79289637fa0ef6d4436ebdfe6905cbd6"
uuid = "e2d170a0-9d28-54be-80f0-106bbe20a464"
version = "1.0.0"

[[DataValues]]
deps = ["DataValueInterfaces", "Dates"]
git-tree-sha1 = "d88a19299eba280a6d062e135a43f00323ae70bf"
uuid = "e7dc6d0d-1eca-5fa6-8ad6-5aecde8b7ea5"
version = "0.4.13"

[[Dates]]
deps = ["Printf"]
uuid = "ade2ca70-3891-5945-98fb-dc099432e06a"

[[DelimitedFiles]]
deps = ["Mmap"]
uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab"

[[Distances]]
deps = ["LinearAlgebra", "Statistics", "StatsAPI"]
git-tree-sha1 = "9f46deb4d4ee4494ffb5a40a27a2aced67bdd838"
uuid = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7"
version = "0.10.4"

[[Distributed]]
deps = ["Random", "Serialization", "Sockets"]
uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"

[[Distributions]]
deps = ["ChainRulesCore", "FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns"]
git-tree-sha1 = "a9b99024b57d12fb19892d3f2230856f6d9671a4"
uuid = "31c24e10-a181-5473-b8eb-7969acd0382f"
version = "0.25.17"

[[DocStringExtensions]]
deps = ["LibGit2"]
git-tree-sha1 = "a32185f5428d3986f47c2ab78b1f216d5e6cc96f"
uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
version = "0.8.5"

[[Downloads]]
deps = ["ArgTools", "LibCURL", "NetworkOptions"]
uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"

[[EarCut_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "3f3a2501fa7236e9b911e0f7a588c657e822bb6d"
uuid = "5ae413db-bbd1-5e63-b57d-d24a61df00f5"
version = "2.2.3+0"

[[Expat_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "b3bfd02e98aedfa5cf885665493c5598c350cd2f"
uuid = "2e619515-83b5-522b-bb60-26c02a35a201"
version = "2.2.10+0"

[[FFMPEG]]
deps = ["FFMPEG_jll"]
git-tree-sha1 = "b57e3acbe22f8484b4b5ff66a7499717fe1a9cc8"
uuid = "c87230d0-a227-11e9-1b43-d7ebe4e7570a"
version = "0.4.1"

[[FFMPEG_jll]]
deps = ["Artifacts", "Bzip2_jll", "FreeType2_jll", "FriBidi_jll", "JLLWrappers", "LAME_jll", "Libdl", "Ogg_jll", "OpenSSL_jll", "Opus_jll", "Pkg", "Zlib_jll", "libass_jll", "libfdk_aac_jll", "libvorbis_jll", "x264_jll", "x265_jll"]
git-tree-sha1 = "d8a578692e3077ac998b50c0217dfd67f21d1e5f"
uuid = "b22a6f82-2f65-5046-a5b2-351ab43fb4e5"
version = "4.4.0+0"

[[FFTW]]
deps = ["AbstractFFTs", "FFTW_jll", "LinearAlgebra", "MKL_jll", "Preferences", "Reexport"]
git-tree-sha1 = "463cb335fa22c4ebacfd1faba5fde14edb80d96c"
uuid = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
version = "1.4.5"

[[FFTW_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "c6033cc3892d0ef5bb9cd29b7f2f0331ea5184ea"
uuid = "f5851436-0d7a-5f13-b9de-f02708fd171a"
version = "3.3.10+0"

[[FillArrays]]
deps = ["LinearAlgebra", "Random", "SparseArrays", "Statistics"]
git-tree-sha1 = "29890dfbc427afa59598b8cfcc10034719bd7744"
uuid = "1a297f60-69ca-5386-bcde-b61e274b549b"
version = "0.12.6"

[[FixedPointNumbers]]
deps = ["Statistics"]
git-tree-sha1 = "335bfdceacc84c5cdf16aadc768aa5ddfc5383cc"
uuid = "53c48c17-4a7d-5ca2-90c5-79b7896eea93"
version = "0.8.4"

[[Fontconfig_jll]]
deps = ["Artifacts", "Bzip2_jll", "Expat_jll", "FreeType2_jll", "JLLWrappers", "Libdl", "Libuuid_jll", "Pkg", "Zlib_jll"]
git-tree-sha1 = "21efd19106a55620a188615da6d3d06cd7f6ee03"
uuid = "a3f928ae-7b40-5064-980b-68af3947d34b"
version = "2.13.93+0"

[[Formatting]]
deps = ["Printf"]
git-tree-sha1 = "8339d61043228fdd3eb658d86c926cb282ae72a8"
uuid = "59287772-0a20-5a39-b81b-1366585eb4c0"
version = "0.4.2"

[[FreeType2_jll]]
deps = ["Artifacts", "Bzip2_jll", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
git-tree-sha1 = "87eb71354d8ec1a96d4a7636bd57a7347dde3ef9"
uuid = "d7e528f0-a631-5988-bf34-fe36492bcfd7"
version = "2.10.4+0"

[[FriBidi_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "aa31987c2ba8704e23c6c8ba8a4f769d5d7e4f91"
uuid = "559328eb-81f9-559d-9380-de523a88c83c"
version = "1.0.10+0"

[[GLFW_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libglvnd_jll", "Pkg", "Xorg_libXcursor_jll", "Xorg_libXi_jll", "Xorg_libXinerama_jll", "Xorg_libXrandr_jll"]
git-tree-sha1 = "dba1e8614e98949abfa60480b13653813d8f0157"
uuid = "0656b61e-2033-5cc2-a64a-77c0f6c09b89"
version = "3.3.5+0"

[[GR]]
deps = ["Base64", "DelimitedFiles", "GR_jll", "HTTP", "JSON", "Libdl", "LinearAlgebra", "Pkg", "Printf", "Random", "Serialization", "Sockets", "Test", "UUIDs"]
git-tree-sha1 = "c2178cfbc0a5a552e16d097fae508f2024de61a3"
uuid = "28b8d3ca-fb5f-59d9-8090-bfdbd6d07a71"
version = "0.59.0"

[[GR_jll]]
deps = ["Artifacts", "Bzip2_jll", "Cairo_jll", "FFMPEG_jll", "Fontconfig_jll", "GLFW_jll", "JLLWrappers", "JpegTurbo_jll", "Libdl", "Libtiff_jll", "Pixman_jll", "Pkg", "Qt5Base_jll", "Zlib_jll", "libpng_jll"]
git-tree-sha1 = "ef49a187604f865f4708c90e3f431890724e9012"
uuid = "d2c73de3-f751-5644-a686-071e5b155ba9"
version = "0.59.0+0"

[[GeometryBasics]]
deps = ["EarCut_jll", "IterTools", "LinearAlgebra", "StaticArrays", "StructArrays", "Tables"]
git-tree-sha1 = "58bcdf5ebc057b085e58d95c138725628dd7453c"
uuid = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
version = "0.4.1"

[[Gettext_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "XML2_jll"]
git-tree-sha1 = "9b02998aba7bf074d14de89f9d37ca24a1a0b046"
uuid = "78b55507-aeef-58d4-861c-77aaff3498b1"
version = "0.21.0+0"

[[Glib_jll]]
deps = ["Artifacts", "Gettext_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Libiconv_jll", "Libmount_jll", "PCRE_jll", "Pkg", "Zlib_jll"]
git-tree-sha1 = "7bf67e9a481712b3dbe9cb3dac852dc4b1162e02"
uuid = "7746bdde-850d-59dc-9ae8-88ece973131d"
version = "2.68.3+0"

[[Graphite2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "344bf40dcab1073aca04aa0df4fb092f920e4011"
uuid = "3b182d85-2403-5c21-9c21-1e1f0cc25472"
version = "1.3.14+0"

[[Grisu]]
git-tree-sha1 = "53bb909d1151e57e2484c3d1b53e19552b887fb2"
uuid = "42e2da0e-8278-4e71-bc24-59509adca0fe"
version = "1.0.2"

[[HTTP]]
deps = ["Base64", "Dates", "IniFile", "Logging", "MbedTLS", "NetworkOptions", "Sockets", "URIs"]
git-tree-sha1 = "14eece7a3308b4d8be910e265c724a6ba51a9798"
uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"
version = "0.9.16"

[[HarfBuzz_jll]]
deps = ["Artifacts", "Cairo_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "Graphite2_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Pkg"]
git-tree-sha1 = "8a954fed8ac097d5be04921d595f741115c1b2ad"
uuid = "2e76f6c2-a576-52d4-95c1-20adfe4de566"
version = "2.8.1+0"

[[HypertextLiteral]]
git-tree-sha1 = "72053798e1be56026b81d4e2682dbe58922e5ec9"
uuid = "ac1192a8-f4b3-4bfe-ba22-af5b92cd3ab2"
version = "0.9.0"

[[IOCapture]]
deps = ["Logging", "Random"]
git-tree-sha1 = "f7be53659ab06ddc986428d3a9dcc95f6fa6705a"
uuid = "b5f81e59-6552-4d32-b1f0-c071b021bf89"
version = "0.2.2"

[[IniFile]]
deps = ["Test"]
git-tree-sha1 = "098e4d2c533924c921f9f9847274f2ad89e018b8"
uuid = "83e8ac13-25f8-5344-8a64-a9f2b223428f"
version = "0.5.0"

[[IntelOpenMP_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "d979e54b71da82f3a65b62553da4fc3d18c9004c"
uuid = "1d5cc7b8-4909-519e-a0f8-d0f5ad9712d0"
version = "2018.0.3+2"

[[InteractiveUtils]]
deps = ["Markdown"]
uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240"

[[Interpolations]]
deps = ["AxisAlgorithms", "ChainRulesCore", "LinearAlgebra", "OffsetArrays", "Random", "Ratios", "Requires", "SharedArrays", "SparseArrays", "StaticArrays", "WoodburyMatrices"]
git-tree-sha1 = "61aa005707ea2cebf47c8d780da8dc9bc4e0c512"
uuid = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59"
version = "0.13.4"

[[IrrationalConstants]]
git-tree-sha1 = "f76424439413893a832026ca355fe273e93bce94"
uuid = "92d709cd-6900-40b7-9082-c6be49f344b6"
version = "0.1.0"

[[IterTools]]
git-tree-sha1 = "05110a2ab1fc5f932622ffea2a003221f4782c18"
uuid = "c8e1da08-722c-5040-9ed9-7db0dc04731e"
version = "1.3.0"

[[IteratorInterfaceExtensions]]
git-tree-sha1 = "a3f24677c21f5bbe9d2a714f95dcd58337fb2856"
uuid = "82899510-4779-5014-852e-03e436cf321d"
version = "1.0.0"

[[JLLWrappers]]
deps = ["Preferences"]
git-tree-sha1 = "642a199af8b68253517b80bd3bfd17eb4e84df6e"
uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210"
version = "1.3.0"

[[JSON]]
deps = ["Dates", "Mmap", "Parsers", "Unicode"]
git-tree-sha1 = "8076680b162ada2a031f707ac7b4953e30667a37"
uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
version = "0.21.2"

[[JpegTurbo_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "d735490ac75c5cb9f1b00d8b5509c11984dc6943"
uuid = "aacddb02-875f-59d6-b918-886e6ef4fbf8"
version = "2.1.0+0"

[[KernelDensity]]
deps = ["Distributions", "DocStringExtensions", "FFTW", "Interpolations", "StatsBase"]
git-tree-sha1 = "591e8dc09ad18386189610acafb970032c519707"
uuid = "5ab0869b-81aa-558d-bb23-cbf5423bbe9b"
version = "0.6.3"

[[LAME_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "f6250b16881adf048549549fba48b1161acdac8c"
uuid = "c1c5ebd0-6772-5130-a774-d5fcae4a789d"
version = "3.100.1+0"

[[LZO_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "e5b909bcf985c5e2605737d2ce278ed791b89be6"
uuid = "dd4b983a-f0e5-5f8d-a1b7-129d4a5fb1ac"
version = "2.10.1+0"

[[LaTeXStrings]]
git-tree-sha1 = "c7f1c695e06c01b95a67f0cd1d34994f3e7db104"
uuid = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
version = "1.2.1"

[[Latexify]]
deps = ["Formatting", "InteractiveUtils", "LaTeXStrings", "MacroTools", "Markdown", "Printf", "Requires"]
git-tree-sha1 = "a4b12a1bd2ebade87891ab7e36fdbce582301a92"
uuid = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
version = "0.15.6"

[[LazyArtifacts]]
deps = ["Artifacts", "Pkg"]
uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3"

[[LibCURL]]
deps = ["LibCURL_jll", "MozillaCACerts_jll"]
uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21"

[[LibCURL_jll]]
deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"]
uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0"

[[LibGit2]]
deps = ["Base64", "NetworkOptions", "Printf", "SHA"]
uuid = "76f85450-5226-5b5a-8eaa-529ad045b433"

[[LibSSH2_jll]]
deps = ["Artifacts", "Libdl", "MbedTLS_jll"]
uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8"

[[Libdl]]
uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb"

[[Libffi_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "761a393aeccd6aa92ec3515e428c26bf99575b3b"
uuid = "e9f186c6-92d2-5b65-8a66-fee21dc1b490"
version = "3.2.2+0"

[[Libgcrypt_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libgpg_error_jll", "Pkg"]
git-tree-sha1 = "64613c82a59c120435c067c2b809fc61cf5166ae"
uuid = "d4300ac3-e22c-5743-9152-c294e39db1e4"
version = "1.8.7+0"

[[Libglvnd_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll", "Xorg_libXext_jll"]
git-tree-sha1 = "7739f837d6447403596a75d19ed01fd08d6f56bf"
uuid = "7e76a0d4-f3c7-5321-8279-8d96eeed0f29"
version = "1.3.0+3"

[[Libgpg_error_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "c333716e46366857753e273ce6a69ee0945a6db9"
uuid = "7add5ba3-2f88-524e-9cd5-f83b8a55f7b8"
version = "1.42.0+0"

[[Libiconv_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "42b62845d70a619f063a7da093d995ec8e15e778"
uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531"
version = "1.16.1+1"

[[Libmount_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "9c30530bf0effd46e15e0fdcf2b8636e78cbbd73"
uuid = "4b2f31a3-9ecc-558c-b454-b3730dcb73e9"
version = "2.35.0+0"

[[Libtiff_jll]]
deps = ["Artifacts", "JLLWrappers", "JpegTurbo_jll", "Libdl", "Pkg", "Zlib_jll", "Zstd_jll"]
git-tree-sha1 = "340e257aada13f95f98ee352d316c3bed37c8ab9"
uuid = "89763e89-9b03-5906-acba-b20f662cd828"
version = "4.3.0+0"

[[Libuuid_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "7f3efec06033682db852f8b3bc3c1d2b0a0ab066"
uuid = "38a345b3-de98-5d2b-a5d3-14cd9215e700"
version = "2.36.0+0"

[[LinearAlgebra]]
deps = ["Libdl"]
uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

[[LogExpFunctions]]
deps = ["ChainRulesCore", "DocStringExtensions", "IrrationalConstants", "LinearAlgebra"]
git-tree-sha1 = "34dc30f868e368f8a17b728a1238f3fcda43931a"
uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
version = "0.3.3"

[[Logging]]
uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"

[[MKL_jll]]
deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "Pkg"]
git-tree-sha1 = "5455aef09b40e5020e1520f551fa3135040d4ed0"
uuid = "856f044c-d86e-5d09-b602-aeab76dc8ba7"
version = "2021.1.1+2"

[[MacroTools]]
deps = ["Markdown", "Random"]
git-tree-sha1 = "5a5bc6bf062f0f95e62d0fe0a2d99699fed82dd9"
uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
version = "0.5.8"

[[Markdown]]
deps = ["Base64"]
uuid = "d6f4376e-aef5-505a-96c1-9c027394607a"

[[MbedTLS]]
deps = ["Dates", "MbedTLS_jll", "Random", "Sockets"]
git-tree-sha1 = "1c38e51c3d08ef2278062ebceade0e46cefc96fe"
uuid = "739be429-bea8-5141-9913-cc70e7f3736d"
version = "1.0.3"

[[MbedTLS_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1"

[[Measures]]
git-tree-sha1 = "e498ddeee6f9fdb4551ce855a46f54dbd900245f"
uuid = "442fdcdd-2543-5da2-b0f3-8c86c306513e"
version = "0.3.1"

[[Missings]]
deps = ["DataAPI"]
git-tree-sha1 = "bf210ce90b6c9eed32d25dbcae1ebc565df2687f"
uuid = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28"
version = "1.0.2"

[[Mmap]]
uuid = "a63ad114-7e13-5084-954f-fe012c677804"

[[MozillaCACerts_jll]]
uuid = "14a3606d-f60d-562e-9121-12d972cd8159"

[[MultivariateStats]]
deps = ["Arpack", "LinearAlgebra", "SparseArrays", "Statistics", "StatsBase"]
git-tree-sha1 = "8d958ff1854b166003238fe191ec34b9d592860a"
uuid = "6f286f6a-111f-5878-ab1e-185364afe411"
version = "0.8.0"

[[NaNMath]]
git-tree-sha1 = "bfe47e760d60b82b66b61d2d44128b62e3a369fb"
uuid = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3"
version = "0.3.5"

[[NearestNeighbors]]
deps = ["Distances", "StaticArrays"]
git-tree-sha1 = "16baacfdc8758bc374882566c9187e785e85c2f0"
uuid = "b8a86587-4115-5ab1-83bc-aa920d37bbce"
version = "0.4.9"

[[NetworkOptions]]
uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908"

[[Observables]]
git-tree-sha1 = "fe29afdef3d0c4a8286128d4e45cc50621b1e43d"
uuid = "510215fc-4207-5dde-b226-833fc4488ee2"
version = "0.4.0"

[[OffsetArrays]]
deps = ["Adapt"]
git-tree-sha1 = "c0e9e582987d36d5a61e650e6e543b9e44d9914b"
uuid = "6fe1bfb0-de20-5000-8ca7-80f57d26f881"
version = "1.10.7"

[[Ogg_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "7937eda4681660b4d6aeeecc2f7e1c81c8ee4e2f"
uuid = "e7412a2a-1a6e-54c0-be00-318e2571c051"
version = "1.3.5+0"

[[OpenBLAS_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"]
uuid = "4536629a-c528-5b80-bd46-f80d51c5b363"

[[OpenLibm_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "05823500-19ac-5b8b-9628-191a04bc5112"

[[OpenSSL_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "15003dcb7d8db3c6c857fda14891a539a8f2705a"
uuid = "458c3c95-2e84-50aa-8efc-19380b2a3a95"
version = "1.1.10+0"

[[OpenSpecFun_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1"
uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e"
version = "0.5.5+0"

[[Opus_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "51a08fb14ec28da2ec7a927c4337e4332c2a4720"
uuid = "91d4177d-7536-5919-b921-800302f37372"
version = "1.3.2+0"

[[OrderedCollections]]
git-tree-sha1 = "85f8e6578bf1f9ee0d11e7bb1b1456435479d47c"
uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
version = "1.4.1"

[[PCRE_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "b2a7af664e098055a7529ad1a900ded962bca488"
uuid = "2f80f16e-611a-54ab-bc61-aa92de5b98fc"
version = "8.44.0+0"

[[PDMats]]
deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"]
git-tree-sha1 = "4dd403333bcf0909341cfe57ec115152f937d7d8"
uuid = "90014a1f-27ba-587c-ab20-58faa44d9150"
version = "0.11.1"

[[Parsers]]
deps = ["Dates"]
git-tree-sha1 = "a8709b968a1ea6abc2dc1967cb1db6ac9a00dfb6"
uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
version = "2.0.5"

[[Pixman_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "b4f5d02549a10e20780a24fce72bea96b6329e29"
uuid = "30392449-352a-5448-841d-b1acce4e97dc"
version = "0.40.1+0"

[[Pkg]]
deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"]
uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"

[[PlotThemes]]
deps = ["PlotUtils", "Requires", "Statistics"]
git-tree-sha1 = "a3a964ce9dc7898193536002a6dd892b1b5a6f1d"
uuid = "ccf2f8ad-2431-5c83-bf29-c5338b663b6a"
version = "2.0.1"

[[PlotUtils]]
deps = ["ColorSchemes", "Colors", "Dates", "Printf", "Random", "Reexport", "Statistics"]
git-tree-sha1 = "2537ed3c0ed5e03896927187f5f2ee6a4ab342db"
uuid = "995b91a9-d308-5afd-9ec6-746e21dbc043"
version = "1.0.14"

[[Plots]]
deps = ["Base64", "Contour", "Dates", "Downloads", "FFMPEG", "FixedPointNumbers", "GR", "GeometryBasics", "JSON", "Latexify", "LinearAlgebra", "Measures", "NaNMath", "PlotThemes", "PlotUtils", "Printf", "REPL", "Random", "RecipesBase", "RecipesPipeline", "Reexport", "Requires", "Scratch", "Showoff", "SparseArrays", "Statistics", "StatsBase", "UUIDs"]
git-tree-sha1 = "cfbd033def161db9494f86c5d18fbf874e09e514"
uuid = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
version = "1.22.3"

[[PlutoUI]]
deps = ["Base64", "Dates", "HypertextLiteral", "IOCapture", "InteractiveUtils", "JSON", "Logging", "Markdown", "Random", "Reexport", "UUIDs"]
git-tree-sha1 = "d1fb76655a95bf6ea4348d7197b22e889a4375f4"
uuid = "7f904dfe-b85e-4ff6-b463-dae2292396a8"
version = "0.7.14"

[[Preferences]]
deps = ["TOML"]
git-tree-sha1 = "00cfd92944ca9c760982747e9a1d0d5d86ab1e5a"
uuid = "21216c6a-2e73-6563-6e65-726566657250"
version = "1.2.2"

[[Printf]]
deps = ["Unicode"]
uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7"

[[Qt5Base_jll]]
deps = ["Artifacts", "CompilerSupportLibraries_jll", "Fontconfig_jll", "Glib_jll", "JLLWrappers", "Libdl", "Libglvnd_jll", "OpenSSL_jll", "Pkg", "Xorg_libXext_jll", "Xorg_libxcb_jll", "Xorg_xcb_util_image_jll", "Xorg_xcb_util_keysyms_jll", "Xorg_xcb_util_renderutil_jll", "Xorg_xcb_util_wm_jll", "Zlib_jll", "xkbcommon_jll"]
git-tree-sha1 = "ad368663a5e20dbb8d6dc2fddeefe4dae0781ae8"
uuid = "ea2cea3b-5b76-57ae-a6ef-0a8af62496e1"
version = "5.15.3+0"

[[QuadGK]]
deps = ["DataStructures", "LinearAlgebra"]
git-tree-sha1 = "78aadffb3efd2155af139781b8a8df1ef279ea39"
uuid = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
version = "2.4.2"

[[REPL]]
deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"]
uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb"

[[Random]]
deps = ["Serialization"]
uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

[[Ratios]]
deps = ["Requires"]
git-tree-sha1 = "01d341f502250e81f6fec0afe662aa861392a3aa"
uuid = "c84ed2f1-dad5-54f0-aa8e-dbefe2724439"
version = "0.4.2"

[[RecipesBase]]
git-tree-sha1 = "44a75aa7a527910ee3d1751d1f0e4148698add9e"
uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
version = "1.1.2"

[[RecipesPipeline]]
deps = ["Dates", "NaNMath", "PlotUtils", "RecipesBase"]
git-tree-sha1 = "7ad0dfa8d03b7bcf8c597f59f5292801730c55b8"
uuid = "01d81517-befc-4cb6-b9ec-a95719d0359c"
version = "0.4.1"

[[Reexport]]
git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b"
uuid = "189a3867-3050-52da-a836-e630ba90ab69"
version = "1.2.2"

[[Requires]]
deps = ["UUIDs"]
git-tree-sha1 = "4036a3bd08ac7e968e27c203d45f5fff15020621"
uuid = "ae029012-a4dd-5104-9daa-d747884805df"
version = "1.1.3"

[[Rmath]]
deps = ["Random", "Rmath_jll"]
git-tree-sha1 = "bf3188feca147ce108c76ad82c2792c57abe7b1f"
uuid = "79098fc4-a85e-5d69-aa6a-4863f24498fa"
version = "0.7.0"

[[Rmath_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "68db32dff12bb6127bac73c209881191bf0efbb7"
uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f"
version = "0.3.0+0"

[[SHA]]
uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"

[[Scratch]]
deps = ["Dates"]
git-tree-sha1 = "0b4b7f1393cff97c33891da2a0bf69c6ed241fda"
uuid = "6c6a2e73-6563-6170-7368-637461726353"
version = "1.1.0"

[[SentinelArrays]]
deps = ["Dates", "Random"]
git-tree-sha1 = "54f37736d8934a12a200edea2f9206b03bdf3159"
uuid = "91c51154-3ec4-41a3-a24f-3f23e20d615c"
version = "1.3.7"

[[Serialization]]
uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"

[[SharedArrays]]
deps = ["Distributed", "Mmap", "Random", "Serialization"]
uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383"

[[Showoff]]
deps = ["Dates", "Grisu"]
git-tree-sha1 = "91eddf657aca81df9ae6ceb20b959ae5653ad1de"
uuid = "992d4aef-0814-514b-bc4d-f2e9a6c4116f"
version = "1.0.3"

[[Sockets]]
uuid = "6462fe0b-24de-5631-8697-dd941f90decc"

[[SortingAlgorithms]]
deps = ["DataStructures"]
git-tree-sha1 = "b3363d7460f7d098ca0912c69b082f75625d7508"
uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c"
version = "1.0.1"

[[SparseArrays]]
deps = ["LinearAlgebra", "Random"]
uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[[SpecialFunctions]]
deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"]
git-tree-sha1 = "793793f1df98e3d7d554b65a107e9c9a6399a6ed"
uuid = "276daf66-3868-5448-9aa4-cd146d93841b"
version = "1.7.0"

[[StaticArrays]]
deps = ["LinearAlgebra", "Random", "Statistics"]
git-tree-sha1 = "3c76dde64d03699e074ac02eb2e8ba8254d428da"
uuid = "90137ffa-7385-5640-81b9-e52037218182"
version = "1.2.13"

[[Statistics]]
deps = ["LinearAlgebra", "SparseArrays"]
uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[[StatsAPI]]
git-tree-sha1 = "1958272568dc176a1d881acb797beb909c785510"
uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0"
version = "1.0.0"

[[StatsBase]]
deps = ["DataAPI", "DataStructures", "LinearAlgebra", "Missings", "Printf", "Random", "SortingAlgorithms", "SparseArrays", "Statistics", "StatsAPI"]
git-tree-sha1 = "8cbbc098554648c84f79a463c9ff0fd277144b6c"
uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
version = "0.33.10"

[[StatsFuns]]
deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"]
git-tree-sha1 = "95072ef1a22b057b1e80f73c2a89ad238ae4cfff"
uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c"
version = "0.9.12"

[[StatsPlots]]
deps = ["Clustering", "DataStructures", "DataValues", "Distributions", "Interpolations", "KernelDensity", "LinearAlgebra", "MultivariateStats", "Observables", "Plots", "RecipesBase", "RecipesPipeline", "Reexport", "StatsBase", "TableOperations", "Tables", "Widgets"]
git-tree-sha1 = "eb007bb78d8a46ab98cd14188e3cec139a4476cf"
uuid = "f3b207a7-027a-5e70-b257-86293d7955fd"
version = "0.14.28"

[[StructArrays]]
deps = ["Adapt", "DataAPI", "StaticArrays", "Tables"]
git-tree-sha1 = "2ce41e0d042c60ecd131e9fb7154a3bfadbf50d3"
uuid = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
version = "0.6.3"

[[SuiteSparse]]
deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"]
uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9"

[[TOML]]
deps = ["Dates"]
uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76"

[[TableOperations]]
deps = ["SentinelArrays", "Tables", "Test"]
git-tree-sha1 = "019acfd5a4a6c5f0f38de69f2ff7ed527f1881da"
uuid = "ab02a1b2-a7df-11e8-156e-fb1833f50b87"
version = "1.1.0"

[[TableTraits]]
deps = ["IteratorInterfaceExtensions"]
git-tree-sha1 = "c06b2f539df1c6efa794486abfb6ed2022561a39"
uuid = "3783bdb8-4a98-5b6b-af9a-565f29a5fe9c"
version = "1.0.1"

[[Tables]]
deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "TableTraits", "Test"]
git-tree-sha1 = "1162ce4a6c4b7e31e0e6b14486a6986951c73be9"
uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
version = "1.5.2"

[[Tar]]
deps = ["ArgTools", "SHA"]
uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e"

[[Test]]
deps = ["InteractiveUtils", "Logging", "Random", "Serialization"]
uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

[[URIs]]
git-tree-sha1 = "97bbe755a53fe859669cd907f2d96aee8d2c1355"
uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4"
version = "1.3.0"

[[UUIDs]]
deps = ["Random", "SHA"]
uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"

[[Unicode]]
uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5"

[[Wayland_jll]]
deps = ["Artifacts", "Expat_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Pkg", "XML2_jll"]
git-tree-sha1 = "3e61f0b86f90dacb0bc0e73a0c5a83f6a8636e23"
uuid = "a2964d1f-97da-50d4-b82a-358c7fce9d89"
version = "1.19.0+0"

[[Wayland_protocols_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Wayland_jll"]
git-tree-sha1 = "2839f1c1296940218e35df0bbb220f2a79686670"
uuid = "2381bf8a-dfd0-557d-9999-79630e7b1b91"
version = "1.18.0+4"

[[Widgets]]
deps = ["Colors", "Dates", "Observables", "OrderedCollections"]
git-tree-sha1 = "80661f59d28714632132c73779f8becc19a113f2"
uuid = "cc8bc4a8-27d6-5769-a93b-9d913e69aa62"
version = "0.6.4"

[[WoodburyMatrices]]
deps = ["LinearAlgebra", "SparseArrays"]
git-tree-sha1 = "59e2ad8fd1591ea019a5259bd012d7aee15f995c"
uuid = "efce3f68-66dc-5838-9240-27a6d6f5f9b6"
version = "0.5.3"

[[XML2_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"]
git-tree-sha1 = "1acf5bdf07aa0907e0a37d3718bb88d4b687b74a"
uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a"
version = "2.9.12+0"

[[XSLT_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libgcrypt_jll", "Libgpg_error_jll", "Libiconv_jll", "Pkg", "XML2_jll", "Zlib_jll"]
git-tree-sha1 = "91844873c4085240b95e795f692c4cec4d805f8a"
uuid = "aed1982a-8fda-507f-9586-7b0439959a61"
version = "1.1.34+0"

[[Xorg_libX11_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libxcb_jll", "Xorg_xtrans_jll"]
git-tree-sha1 = "5be649d550f3f4b95308bf0183b82e2582876527"
uuid = "4f6342f7-b3d2-589e-9d20-edeb45f2b2bc"
version = "1.6.9+4"

[[Xorg_libXau_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "4e490d5c960c314f33885790ed410ff3a94ce67e"
uuid = "0c0b7dd1-d40b-584c-a123-a41640f87eec"
version = "1.0.9+4"

[[Xorg_libXcursor_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXfixes_jll", "Xorg_libXrender_jll"]
git-tree-sha1 = "12e0eb3bc634fa2080c1c37fccf56f7c22989afd"
uuid = "935fb764-8cf2-53bf-bb30-45bb1f8bf724"
version = "1.2.0+4"

[[Xorg_libXdmcp_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "4fe47bd2247248125c428978740e18a681372dd4"
uuid = "a3789734-cfe1-5b06-b2d0-1dd0d9d62d05"
version = "1.1.3+4"

[[Xorg_libXext_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "b7c0aa8c376b31e4852b360222848637f481f8c3"
uuid = "1082639a-0dae-5f34-9b06-72781eeb8cb3"
version = "1.3.4+4"

[[Xorg_libXfixes_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "0e0dc7431e7a0587559f9294aeec269471c991a4"
uuid = "d091e8ba-531a-589c-9de9-94069b037ed8"
version = "5.0.3+4"

[[Xorg_libXi_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXext_jll", "Xorg_libXfixes_jll"]
git-tree-sha1 = "89b52bc2160aadc84d707093930ef0bffa641246"
uuid = "a51aa0fd-4e3c-5386-b890-e753decda492"
version = "1.7.10+4"

[[Xorg_libXinerama_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXext_jll"]
git-tree-sha1 = "26be8b1c342929259317d8b9f7b53bf2bb73b123"
uuid = "d1454406-59df-5ea1-beac-c340f2130bc3"
version = "1.1.4+4"

[[Xorg_libXrandr_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libXext_jll", "Xorg_libXrender_jll"]
git-tree-sha1 = "34cea83cb726fb58f325887bf0612c6b3fb17631"
uuid = "ec84b674-ba8e-5d96-8ba1-2a689ba10484"
version = "1.5.2+4"

[[Xorg_libXrender_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "19560f30fd49f4d4efbe7002a1037f8c43d43b96"
uuid = "ea2f1a96-1ddc-540d-b46f-429655e07cfa"
version = "0.9.10+4"

[[Xorg_libpthread_stubs_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "6783737e45d3c59a4a4c4091f5f88cdcf0908cbb"
uuid = "14d82f49-176c-5ed1-bb49-ad3f5cbd8c74"
version = "0.1.0+3"

[[Xorg_libxcb_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "XSLT_jll", "Xorg_libXau_jll", "Xorg_libXdmcp_jll", "Xorg_libpthread_stubs_jll"]
git-tree-sha1 = "daf17f441228e7a3833846cd048892861cff16d6"
uuid = "c7cfdc94-dc32-55de-ac96-5a1b8d977c5b"
version = "1.13.0+3"

[[Xorg_libxkbfile_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libX11_jll"]
git-tree-sha1 = "926af861744212db0eb001d9e40b5d16292080b2"
uuid = "cc61e674-0454-545c-8b26-ed2c68acab7a"
version = "1.1.0+4"

[[Xorg_xcb_util_image_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "0fab0a40349ba1cba2c1da699243396ff8e94b97"
uuid = "12413925-8142-5f55-bb0e-6d7ca50bb09b"
version = "0.4.0+1"

[[Xorg_xcb_util_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libxcb_jll"]
git-tree-sha1 = "e7fd7b2881fa2eaa72717420894d3938177862d1"
uuid = "2def613f-5ad1-5310-b15b-b15d46f528f5"
version = "0.4.0+1"

[[Xorg_xcb_util_keysyms_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "d1151e2c45a544f32441a567d1690e701ec89b00"
uuid = "975044d2-76e6-5fbe-bf08-97ce7c6574c7"
version = "0.4.0+1"

[[Xorg_xcb_util_renderutil_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "dfd7a8f38d4613b6a575253b3174dd991ca6183e"
uuid = "0d47668e-0667-5a69-a72c-f761630bfb7e"
version = "0.3.9+1"

[[Xorg_xcb_util_wm_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xcb_util_jll"]
git-tree-sha1 = "e78d10aab01a4a154142c5006ed44fd9e8e31b67"
uuid = "c22f9ab0-d5fe-5066-847c-f4bb1cd4e361"
version = "0.4.1+1"

[[Xorg_xkbcomp_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_libxkbfile_jll"]
git-tree-sha1 = "4bcbf660f6c2e714f87e960a171b119d06ee163b"
uuid = "35661453-b289-5fab-8a00-3d9160c6a3a4"
version = "1.4.2+4"

[[Xorg_xkeyboard_config_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Xorg_xkbcomp_jll"]
git-tree-sha1 = "5c8424f8a67c3f2209646d4425f3d415fee5931d"
uuid = "33bec58e-1273-512f-9401-5d533626f822"
version = "2.27.0+4"

[[Xorg_xtrans_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "79c31e7844f6ecf779705fbc12146eb190b7d845"
uuid = "c5fb5394-a638-5e4d-96e5-b29de1b5cf10"
version = "1.4.0+3"

[[Zlib_jll]]
deps = ["Libdl"]
uuid = "83775a58-1f1d-513f-b197-d71354ab007a"

[[Zstd_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "cc4bf3fdde8b7e3e9fa0351bdeedba1cf3b7f6e6"
uuid = "3161d3a3-bdf6-5164-811a-617609db77b4"
version = "1.5.0+0"

[[libass_jll]]
deps = ["Artifacts", "Bzip2_jll", "FreeType2_jll", "FriBidi_jll", "HarfBuzz_jll", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
git-tree-sha1 = "5982a94fcba20f02f42ace44b9894ee2b140fe47"
uuid = "0ac62f75-1d6f-5e53-bd7c-93b484bb37c0"
version = "0.15.1+0"

[[libfdk_aac_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "daacc84a041563f965be61859a36e17c4e4fcd55"
uuid = "f638f0a6-7fb0-5443-88ba-1cc74229b280"
version = "2.0.2+0"

[[libpng_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
git-tree-sha1 = "94d180a6d2b5e55e447e2d27a29ed04fe79eb30c"
uuid = "b53b4c65-9356-5827-b1ea-8c7a1a84506f"
version = "1.6.38+0"

[[libvorbis_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Ogg_jll", "Pkg"]
git-tree-sha1 = "c45f4e40e7aafe9d086379e5578947ec8b95a8fb"
uuid = "f27f6e37-5d2b-51aa-960f-b287f2bc3b7a"
version = "1.3.7+0"

[[nghttp2_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d"

[[p7zip_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0"

[[x264_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "4fea590b89e6ec504593146bf8b988b2c00922b2"
uuid = "1270edf5-f2f9-52d2-97e9-ab00b5d0237a"
version = "2021.5.5+0"

[[x265_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "ee567a171cce03570d77ad3a43e90218e38937a9"
uuid = "dfaa095f-4041-5dcd-9319-2fabd8486b76"
version = "3.5.0+0"

[[xkbcommon_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Wayland_jll", "Wayland_protocols_jll", "Xorg_libxcb_jll", "Xorg_xkeyboard_config_jll"]
git-tree-sha1 = "ece2350174195bb31de1a63bea3a41ae1aa593b6"
uuid = "d8fb68d0-12a3-5cfd-a85a-d49703b185fd"
version = "0.9.1+5"
"""

# ╔═╡ Cell order:
# ╟─a4a73fb0-2213-11ec-02f8-1d2409575d4c
# ╠═b6681ace-5cb1-49cb-9ed7-e749c7a67428
# ╟─4dcf2480-7a37-4236-b51b-a171407bd8a0
# ╠═e044742f-15dd-42e1-9271-43f11c7e2785
# ╟─c2d41bee-4952-429b-b3d8-aea46eefd2c4
# ╟─ef92132a-211d-4e8f-8ab2-12ae2c4a9c49
# ╠═e2b04571-e4e5-49fe-823d-af0404eceab7
# ╟─e2bb3b47-1973-4d40-91f2-a703470d1a6a
# ╟─8f3de7b2-2fa3-4aa3-ba81-243aa4722a3a
# ╟─1a193ab0-10b0-447d-a68a-0fc978c6e804
# ╟─d142a6c3-e11c-44b5-930f-16deae98e74f
# ╠═33fc930d-8693-41f9-8614-dc71870bd56d
# ╟─4d10cc40-213e-4c4e-bea8-ea318aee18f7
# ╟─9149e9e5-b0b9-4c22-8795-d982a3cacf6a
# ╟─af901869-37bf-4992-bfc3-ee8477a37b34
# ╟─361b07cc-6554-40d1-9448-6f4e4cd4bb71
# ╟─2161be7c-571e-46c9-b3c2-0543d5e3296f
# ╟─e82146d1-ea29-40a3-b89f-5d178f4a900a
# ╟─da10276c-1381-4edb-9171-1ec0820a1809
# ╟─b8e49bec-d582-4695-b7ac-6b0ccbca5e08
# ╟─2cccdace-62df-4771-aaa4-1dd4342335e4
# ╟─a673dd91-ccdf-4783-8b24-357a7307d7d3
# ╟─3a34068e-4fed-40c3-89ef-2aabafb5144e
# ╟─9c599a59-96f7-4b2c-b690-13750a14eef8
# ╟─ffc34a30-9dba-40d0-afca-2279f4753661
# ╟─a70c520f-860f-48be-a1a6-e87607298c37
# ╠═8c493e53-bb7b-4110-a732-1cf93ecb754b
# ╟─9908c97f-efb4-4634-aa97-9b4d2f9d7acc
# ╠═d690e1d2-654c-4da4-9255-d3871429d262
# ╟─ee5bfc95-bbb0-4d1b-8a83-77e884171f05
# ╠═b7033638-277e-4610-89ef-3abe33ef35bf
# ╟─00000000-0000-0000-0000-000000000001
# ╟─00000000-0000-0000-0000-000000000002

max_stars_repo_name: TomRottier/statistical-rethinking

max_issues_repo_name: TomRottier/statistical-rethinking

max_forks_repo_name: TomRottier/statistical-rethinking


text:
#!/usr/bin/env julia

using HDF5
using Cooler

coolName = "MB286.mcool"

startCoord=1000
endCoord = 10000000
resolution = "resolutions/1000"
#resolution = "resolutions/2000"
resolution = "resolutions/16000"

attrs = h5readattr(coolName,resolution)
binsize = attrs["bin-size"]

binStart = startCoord/binsize
binEnd = endCoord/binsize
matwidth = floor(Int,binEnd-binStart+1)

println("binStart: ",binStart,"\tbinEnd: ",binEnd,"\tmatwidth: ",matwidth)

h5open(coolName) do cin
	onek=cin[resolution]
	pixels = onek["pixels"]
	bin1_id = read(pixels,"bin1_id")
	bin2_id = read(pixels,"bin2_id")
	count = read(pixels,"count")

	println("countlen:",length(count))

	matrix = zeros(matwidth,matwidth)
	println("matrix len:",length(matrix))
	nnv = length(count)
	for binIdx in 1:nnv
#	for binIdx in 1:100
		row = bin1_id[binIdx]+1
		col = bin2_id[binIdx]+1

		#println("binIdx:",binIdx,"\trow:",row,"\tcol:",col,"\tcount:",count[binIdx])

		if (row <= matwidth) && (col <= matwidth)
			matrix[row,col] = count[binIdx]
		end
	end
end

#resolution/2000
#binStart: 0.5	binEnd: 5000.0	matwidth: 5000
#countlen:10,350,414
#matrix len:25,000,000

# resolution/1000
# binStart: 1.0	binEnd: 10000.0	matwidth: 10000
# countlen:13,754,067
# matrix len:100,000,000

ext: jl

max_stars_repo_name: jdurbin/sandbox

max_issues_repo_name: jdurbin/sandbox

max_forks_repo_name: jdurbin/sandbox


text:
# setindex.jl
# test

using LinearMapsAA
using Test: @test

"""
`LinearMapAA_test_setindex(A::LinearMapAM)`
"""
function LinearMapAA_test_setindex(A::LinearMapAM)

    @test all(size(A) .>= (4,4)) # required by tests

    tf1 = [false; trues(size(A,1)-1)]
    tf2 = [false; trues(size(A,2)-2); false]
    ii1 = (3, 2:4, [2,4], :, tf1)
    ii2 = (2, 3:4, [1,4], :, tf2)

    # test all [?,?] combinations with array "X"
    for i2 in ii2
        for i1 in ii1
            B = copy(A)
            X = 2 .+ A[i1,i2].^2 # values must differ from A[ii,jj]
            B[i1,i2] = X
            Am = Matrix(A)
            Bm = Matrix(B)
            Am[i1,i2] = X
            @test isapprox(Am, Bm)
        end
    end

    # test all [?,?] combinations with scalar "s"
    for i2 in ii2
        for i1 in ii1
            B = copy(A)
            s = maximum(abs.(A[:])) + 2
            B[i1,i2] = s
            Am = Matrix(A)
            Bm = Matrix(B)
            if (i1 == Colon() || ndims(i1) > 0) || (i2 == Colon() || ndims(i2) > 0)
                Am[i1,i2] .= s
            else
                Am[i1,i2] = s
            end
            @test isapprox(Am, Bm)
        end
    end

    # others not supported for now
    set1 = (3, ) # [3], [2,4], 2:4, (1:length(A)) .== 2), end # todo
    for s1 in set1
        B = copy(A)
        X = 2 .+ A[s1].^2 # values must differ from A[ii,jj]
        B[s1] = X
        Am = Matrix(A)
        Bm = Matrix(B)
        Am[s1] = X
        @test isapprox(Am, Bm)
    end

    # insanity below here

    # A[:] = s
    B = copy(A)
    B[:] = 5
    Am = Matrix(A)
    Am[:] .= 5
    Bm = Matrix(B)
    @test Bm == Am

    # A[:] = v
    B = copy(A)
    v = 1:length(A)
    B[:] = v
    Am = Matrix(A)
    Am[:] .= v
    Bm = Matrix(B)
    @test Bm == Am

    # A[:,:]
    B = copy(A)
    B[:,:] = 6
    Am = Matrix(A)
    Am[:,:] .= 6
    Bm = Matrix(B)
    @test Bm == Am

    true
end

N = 6; M = 8 # non-square to stress test
forw = x -> [cumsum(x); 0; 0]
back = y -> reverse(cumsum(reverse(y[1:N])))
A = LinearMapAA(forw, back, (M, N))
@test LinearMapAA_test_setindex(A)


text:
using MethodOfParticularSolutions
using Nemo
using Test
using Printf

include("sphericaltriangle.jl")
include("lshape.jl")


text:
###########################
# Helper functions for 3 or more group-by byvec
###########################

export diffif, diffif!, meanreduce

function diffif!(difftoprev, vec, l = length(vec))
    lastvec = vec[1]
    @inbounds for i=2:l
        thisvec = vec[i]
        if lastvec != thisvec
            difftoprev[i] = true
            lastvec = thisvec
        end
    end
    difftoprev
end

function diffif(byveccv)
    l = length(byveccv[1])
    diff2prev = BitVector(undef, l)
    diff2prev .= false
    diff2prev[1] = true
    @inbounds for i in 1:length(byveccv)
        diffif!(diff2prev, byveccv[i], l)
    end
    diff2prev
end

###########################
# meanreduce
###########################
meanreduce(mean_so_far, new_value, i) = mean_so_far*(i-1)/i + new_value/i

#############################################
# Single groups tuple of fgroupreduce
# Single fn
# Single val
# as long as grouptwo! is defined for byveccv then it's fine
#############################################
"""
    fgroupreduce(fn, byvec, valvec, init)

Group by `byvec` and apply `reduce(fn, valvec, init = init)` within each group
of `byvec`
"""
fgroupreduce(fn, byvec::AbstractVector{T}, valvec::AbstractVector{Z},  init) where {T, Z} = begin
    cm = Dict{T, typeof(fn(init, valvec[1]))}()
    for (b, v) in zip(byvec, valvec)
        index = ht_keyindex2!(cm, b)
        if index > 0
            @inbounds cm.vals[index] = fn(cm.vals[index], v)
        else
            @inbounds Base._setindex!(cm, fn(init, v), b, -index)
        end
    end
    cm
end

#############################################
# Multiple groups tuple of fgroupreduce
# Single fn
# Single val
# as long as grouptwo! is defined for byveccv then it's fine
#############################################
function fgroupreduce!(fn::F, byveccv, val::Vector{Z}, v0 = zero(Z)) where {F<:Function, Z}
    l = length(val)
    # index = collect(1:l)
    lb = length(byveccv)

    # # sort the first column and index
    # @time grouptwo!(byveccv[lb], index)
    #
    # # reorganise value
    # @time @inbounds val .= val[index]
    #
    # # reorganises the other columns
    # @time @inbounds for i = lb-1:-1:1
    #     byveccv[i] .= byveccv[i][index]
    # end

    c1l = collect(1:l)

    # sort the result of the columns
    @time @inbounds for i = lb:-1:2
        println("sorting $i th column")

        index = copy(c1l)

        grouptwo!(byveccv[i], index)
        # for j = lb:-1:i+1
        #     println("    organising $j th column")
        #     byveccv[j] .= byveccv[j][index]
        # end
        for j = i-1:-1:1
            println("    organising $j th column")
            byveccv[j] .= byveccv[j][index]
        end
        println("    organising value")
        val .= val[index]
    end

    println("sorting 1st column and organising value")
    grouptwo!(byveccv[1], val)

    @time diff2prev = diffif(byveccv)
    n_uniques = sum(diff2prev)

    upto = UInt(0)
    res = fill(v0, n_uniques)

    resby = (bv[diff2prev] for bv in byveccv)

    i = 0
    @time @inbounds for (vali, dp) in zip(val, diff2prev)
        #increase upto by 1 if it's different to previous value
        if UInt(dp) == 1
            i = 0
        end
        i += 1
        upto += UInt(dp)
        #res[upto] = fn(res[upto], vali, i)
        res[upto] = fn(res[upto], vali)
    end
    @time (resby..., res)
end

if false
    fn = +
    x = rand(1:1_000, 100_000_000)
    y = rand(1:1_000, 100_000_000)
    byveccv = (x, y)
    val = copy(x)
    v0 = zero(Int)

     fgroupreduce((a, b, _) ->  a + b, byveccv, val)

    Random.seed!(0)
    z = rand(1:1_000_000, 100_000_000)
    val = copy(z)
     fgroupreduce((cum, new, i)-> cum*(i-1)/i + new/i , z, val, 0.0)

     fgroupreduce((a, b, _) ->  a + b, z, val)

    byveccv = (z,)
    fgroupreduce((cum, new, i)-> cum*(i-1)/i + new/i , x, x, 0.0)
end

# if false
#     byveccv = (categorical(df[:id1]).refs, categorical(df[:id2]).refs) .|> copy
#
#     hehe = DataFrame(deepcopy(collect(byveccv)))
#     sort!(hehe, cols=[:x1,:x2])
#
#     fn = +
#     val = df[:v1]
#     T = Int
#     v0 = 0
#
#     byveccv = (rand(1:100,10_000_000), rand(1:100, 10_000_000))
#     val = rand(1:5,10_000_000)
#
#     df[:id10] =
#
#     byveccv1 = (categorical(df[:id1]).refs, categorical(df[:id2]).refs, rand(1:100, 10_000_000)) .|> copy
#
#     @time FastGroupBy.fgroupreduce!(+, byveccv1, val, 0)
#
#     @time FastGroupBy.fgroupreduce(+, byveccv, val, 0)
#
#     @time FastGroupBy.fgroupreduce!(+, byveccv, val, 0)
#     @time FastGroupBy.fgroupreduce2!(+, byveccv, val, 0)
#
#     @time aggregate(df[[:id1, :id2, :v1]], [:id1,:id2], sum)
#
#     @code_warntype fgroupreduce!(+, byveccv, val, 0)
#
#     @time fgroupreduce(+, byveccv, val, 0)
#
#     @time index = fsortperm(byveccv[2])
#
#     @time v2 = byveccv[1][index]
#     @time index = fsortperm(v2)
# end

#############################################
# Multiple groups tuple of fgroupreduce
# Multiple fn
# Multiple val
# as long as grouptwo! is defined for byveccv then it's fine
# TODO: finish this
#############################################
# function fgroupreduce!(fn::NTuple{M, Function}, byveccv::NTuple{N, AbstractVector}, val::NTuple{M, AbstractVector} , v0 = ((zero(eltype(vt)) for vt in val)...)) where {N, M}
#     lenval = length(val[1])
#     index = collect(1:lenval)
#     grouptwo!(byveccv[N], index)
#
#     # reorders the value vectors
#     @time for i = 1:M
#         @inbounds val[i] .= val[i][index]
#     end
#
#     @time for i = N-1:-1:1
#         @inbounds byveccv[i] .= byveccv[i][index]
#     end
#
#     @time @inbounds for i = N-1:-1:1
#         index .= collect(1:lenval)
#         grouptwo!(byveccv[i], index)
#         for j = N:-1:i+1
#             byveccv[j] .= byveccv[j][index]
#         end
#         for j = 1:M
#             val[j] .= val[j][index]
#         end
#     end
#
#     # diff2prev = diffif(byveccv)
#     # n_uniques = sum(diff2prev)
#
#     # upto::UInt = 0
#     # res = fill(v0, n_uniques)
#
#     # res[1] = v0
#     # resby = ((bv[diff2prev] for bv in byveccv)...)
#     # @inbounds for (vali, dp) in zip(val, diff2prev)
#     #     # increase upto by 1 if it's different to previous value
#     #     upto += UInt(dp)
#     #     res[upto] = fn(res[upto], vali)
#     # end
#     # (resby..., res)
# end
#
# function fgroupreduce(fn::NTuple{M, Function}, byveccv::NTuple{N, AbstractVector}, val::NTuple{M, AbstractVector}, v0 = ((zero(eltype(vt)) for vt in val)...)) where {M, N}
#     fgroupreduce!(fn, ((copy(bv) for bv in byveccv)...), ((copy(v) for v in val)...), v0)
# end
#
# if false
#     fn = +
#     byveccv = (df[:id1], df[:id2])
#     M = 3
#     N= 2
#     val = (df[:v1], df[:v2], df[:v3])
#
#     @time fgroupreduce((+,+,+), byveccv, val)
# end

########################################
# fgroupreduce for single categorical
########################################
function fgroupreduce!(fn::F, byveccv::CategoricalVector, val::Vector{Z}, v0::T = zero(Z)) where {F<:Function, Z,T}
    l1 = length(byveccv.pool)
    lv = length(val)

    # make a histogram of unique values
    res = fill(v0, l1)
    taken = BitArray(l1)
    taken .= false
    @inbounds for i = 1:lv
        k = byveccv.refs[i]
        res[k] = fn(res[k], val[i])
        taken[k] = true
    end

    num_distinct = sum(taken)

    outbyveccv = copy(@view(byveccv[1:num_distinct]))
    outval = Vector{T}(num_distinct)

    distinct_encountered = 1
    @inbounds for i=1:l1
        if taken[i]
            outbyveccv.refs[distinct_encountered] = i
            outval[distinct_encountered] = res[i]
            distinct_encountered += 1
        end
    end

    (outbyveccv, outval)
end

########################################
# fgroupreduce for DataFrames
########################################

# only one group by symbol
fgroupreduce(fn, df::AbstractDataFrame, byveccv::Symbol, val::Symbol, v0 = zero(eltype(df[val]))) =
    DataFrame(fgroupreduce(fn, df[!, byveccv], df[!, val], v0) |> collect, [byveccv, val])

# multiple group by symbol
fgroupreduce(fn, df::AbstractDataFrame, bysyms::NTuple{N, Symbol}, val::Symbol) where N =
DataFrame(
    fgroupreduce(
        fn,
        ((df[bs] for bs in bysyms)...),
        df[val]
    ) |> collect
, [bysyms..., val])

# if false
#     a = "id".*dec.(1:100, 3);
#     ar = rand(a, 10_000_00);
#     val = rand(10_000_00);
#     using FastGroupBy
#     @time fastby(sum, ar, val);
#
#     accv = ar |> CategoricalVector
#
#     @time fgroupreduce(+, accv, val)
#
#     using FastGroupBy
#     @time fastby(sum, a, val)
#
#
#     fgroupreduce(sum, df, :)
# end


text:
function permutation(n1::S, n2::S) where S <: Int

  nn = n1*n2
  p = zeros(S,nn,nn)

  i = 1
  j = 1
  for r = 1:nn
    if j > n2
      j = 1
      i += 1
    end
    c = i+n1*(j-1)
    p[r,c] = 1

    j += 1
  end

  return p

end


text:
export findβ, microcanonical_ensemble, canonical_ensemble

# stupid implementation of softmax - vals should be regularized beforehand
function _softmax(vals)
    w = exp.(vals)
    w ./= sum(w)
end

"""
    _weights_canonical(β, evals)
Compute the distribution of eigenstates with given energies at temperature 1/β.
P(E) ∝ ∑ᵢexp(-β Eᵢ)

Energy values are assumed to be sorted!
"""
function _weights_canonical(evals, β)
    # assume sorted evals!
    regulator = β < 0 ? evals[end] : evals[1]
    _softmax(@. -β * (evals - regulator))
end

_ΔE(evals, E_0, β) = abs((evals ⋅ _weights_canonical(evals, β)) - E_0)

# TODO: How to optimize best? The function is pretty well behaved...
"""
    findβ(evals, E₀; β_0=0)
Compute the inverse temperature β s. t. the thermal expectation value of the energy matches E₀.
"""
findβ(evals, E_0; β_0 = 0) = Optim.optimize(β -> _ΔE(evals, E_0, β[1]), [float(β_0)], Optim.Newton(linesearch = LineSearches.HagerZhang())).minimizer[1]

"""
    canonical_ensemble(evals, E₀; β₀ = 0)
Compute the weights of the state assuming a canonical ensemble of temperature 1/β. β₀ is taken as starting value for the computation.
A state is assigned the weight ∝ exp(-β Eᵢ)

Energy values are assumed to be sorted!
"""
canonical_ensemble(evals, E_0; β_0 = 0) = _weights_canonical(evals, findβ(evals, E_0; β_0))

"""
    microcanonical_ensemble(evals, E₀[, ΔE])
Compute the weights of the state assuming a microcanonical ensemble. Default size of the window is 0.5% of the spectral width.
A state is in the ensemble if |Eᵢ-E_0| < ΔE.
If the ensemble would be empty, choose the closest state (in energy).

Energy values are assumed to be sorted!
"""
function microcanonical_ensemble(evals, E_0, ΔE)::Vector{Float64}
    micro = @. (evals < E_0+ΔE) & (evals > E_0-ΔE)
    s = sum(micro)
    if s > 0
        # enforce float conversion for type stability
        micro ./ float(s)
    else
        micro = zeros(Float64, length(micro))
        smaller = searchsortedlast(evals, E_0)
        ind = abs(evals[smaller]-E_0) < abs(evals[smaller+1]-E_0) ? smaller : smaller+1
        micro[ind] = 1
        micro
    end
end

microcanonical_ensemble(evals, E_0) = microcanonical_ensemble(evals, E_0, 0.005*(evals[end]-evals[1]))


text:
module PlutoLinks

include("./notebook.jl")
include("./revise.jl")
include("./debounce.jl")

export @use_debounce, @use_file, @use_task, @ingredients, @revise

end


text:
export MessageCreate,
    MessageUpdate,
    MessageDelete,
    MessageDeleteBulk,
    MessageReactionAdd,
    MessageReactionRemove,
    MessageReactionRemoveAll

"""
Sent when a [`Message`](@ref) is sent.
"""
struct MessageCreate <: AbstractEvent
    message::Message
end
@boilerplate MessageCreate :docs
MessageCreate(d::Dict{String, Any}) = MessageCreate(Message(d))

"""
Sent when a [`Message`](@ref) is updated.
"""
struct MessageUpdate <: AbstractEvent
    message::Message
end
@boilerplate MessageUpdate :docs
MessageUpdate(d::Dict{String, Any}) = MessageUpdate(Message(d))

"""
Sent when a [`Message`](@ref) is deleted.
"""
struct MessageDelete <: AbstractEvent
    id::Snowflake
    channel_id::Snowflake
    guild_id::Union{Snowflake, Missing}
end
@boilerplate MessageDelete :dict :docs

"""
Sent when multiple [`Message`](@ref)s are deleted in bulk.
"""
struct MessageDeleteBulk <: AbstractEvent
    ids::Vector{Snowflake}
    channel_id::Snowflake
    guild_id::Union{Snowflake, Missing}
end
@boilerplate MessageDeleteBulk :dict :docs

"""
Sent when a [`Reaction`](@ref) is added to a [`Message`](@ref).
"""
struct MessageReactionAdd <: AbstractEvent
    user_id::Snowflake
    channel_id::Snowflake
    message_id::Snowflake
    guild_id::Union{Snowflake, Missing}
    emoji::Emoji
end
@boilerplate MessageReactionAdd :dict :docs

"""
Sent when a [`Reaction`](@ref) is removed from a [`Message`](@ref).
"""
struct MessageReactionRemove <: AbstractEvent
    user_id::Snowflake
    channel_id::Snowflake
    message_id::Snowflake
    guild_id::Union{Snowflake, Missing}
    emoji::Emoji
end
@boilerplate MessageReactionRemove :dict :docs

"""
Sent when all [`Reaction`](@ref)s are removed from a [`Message`](@ref).
"""
struct MessageReactionRemoveAll <: AbstractEvent
    channel_id::Snowflake
    message_id::Snowflake
    guild_id::Union{Snowflake, Missing}
end
@boilerplate MessageReactionRemoveAll :dict :docs


text:
@json struct Committer
    name::String
    date::DateTime
    email::String
end

@json struct CommitTree
    sha::String
    url::String
end

@json struct CommitVerification
    verified::Bool
    reason::String
    signature::String
    payload::String
end

@json struct CommitInfo
    url::String
    author::Committer
    committer::Committer
    message::String
    tree::CommitTree
    comment_count::Int
    verification::CommitVerification
end

@json struct CommitFile
    filename::String
    additions::Int
    deletions::Int
    changes::Int
    status::String
    raw_url::String
    blob_url::String
    patch::String
end

@json struct CommitStats
    additions::Int
    deletions::Int
    total::Int
end

@json struct Commit
    url::String
    sha::String
    node_id::String
    html_url::String
    comments_url::String
    commit::CommitInfo
    author::User
    committer::User
    parents::Vector{CommitTree}
    stats::CommitStats
    files::Vector{CommitFile}
end

endpoint(::GitHubAPI, ::typeof(get_commit), owner::AStr, repo::AStr, ref::AStr) =
    Endpoint(:GET, "/repos/$owner/$repo/commits/$ref")
into(::GitHubAPI, ::typeof(get_commit)) = Commit


text:
"""
    laplacian(n::Int, h::Real)
    laplacian(x::AbstractRange{T}) where {T}

Returns ``n \\times n`` differentiation matrix ``\\mathrm{D^{(2)}}`` using sinc interpolants.
See Weideman, J.A., and Reddy, S.C. (2000). A MATLAB differentiation matrix suite.
ACM Transactions on Mathematical Software (TOMS), 26(4), 465-519.
Equation (20) on page 485.

On input: `n` - size of a uniform grid with step size `h`, or `x` - `AbstractRange` object.

# Example
Second derivative of the gaussian function
```math
    \\partial^2 \\exp(-x^2)/\\partial x^2 = 2(2x^2 - 1)\\exp(-x^2)
```
```jldoctest
julia> x = -6:0.01:6
-6.0:0.01:6.0

julia> d²f_exact = 2(2x.^2 .- 1).*exp.(-x.^2);

julia> d²f_approx = laplacian(length(x), step(x)) * exp.(-x.^2);

julia> isapprox(d²f_exact, d²f_approx, atol = 1e-10)
true
```
"""
function laplacian(n, h)
  Δ = zeros(n, n)
  Δ[diagind(Δ, 0)] .= -1/3*π^2 / h^2
  for i=2:n
    @inbounds Δ[diagind(Δ, i-1)] =
      Δ[diagind(Δ, 1-i)] .= 2*(-1)^i / (i-1)^2 / h^2
  end
  return Δ
end
laplacian(x::AbstractRange) = laplacian(length(x), step(x))

function FGH(x::AbstractRange)
  N, Δx = length(x), step(x)
  iseven(N) && throw(DomainError(N, "number of grid points must be odd"))
  K = zeros(N, N)
  for d=0:N-1
    # See Eq. 26 // C. Clay Marston et al. The Fourier grid Hamiltonian method
    # for bound state eigenvalues and eigenfunctions. JCP 91, 3571 (1989)
    # d ≝ i - j
    K[diagind(K, d)] = K[diagind(K, -d)] .=
    -2/N*(2π/N)^2 / Δx^2 * sum(l -> cos(2π*l*d/N)*l^2, 1:(N >> 1))
  end
  return K
end

function sincdiff(x::AbstractRange)
  n, dx = length(x), step(x)
  ∇ = zeros(n, n)
  ∇[diagind(∇, 0)] .= 0
  for i = 2:n
      @inbounds ∇[diagind(∇, i - 1)] .= (-1)^i / (i - 1) / dx
      @inbounds ∇[diagind(∇, 1 - i)] .= (-1)^(i - 1) / (i - 1) / dx
  end
  return ∇
end

function ∂r(f, x, r = exp.(x), D = sincdiff(x) - I, rₘᵢₙ = 5e-5)
  df = D * (f .* r)
  @. df /= r^2
  extrapolate_left!(df, r, rₘᵢₙ)
  return df
end

#=
 ρ = similar(r)
∇ρ = similar(r)

@. ρ = (2/27*r^2 - 2/3*r + 1)^2 * exp(-2/3*r) / 27π
@. ∇ρ = -(8*r^4 - 192*r^3 + 1512*r^2 - 4536r + 4374)*exp(-2/3*r) / 59049π

norm(∂r(ρ, x) .- ∇ρ)
=#


text:
export MultiAssayExperiment
import SummarizedExperiments
import DataStructures
import DataFrames

function harvest_all_colname(experiments::DataStructures.OrderedDict{String,SummarizedExperiments.SummarizedExperiment})
    all_names = DataStructures.OrderedSet{String}()
    for (key, val) in experiments
        cd = SummarizedExperiments.coldata(val)
        curnames = cd[!,"name"]

        if !isa(curnames, AbstractVector{<:AbstractString})
            throw(ErrorException("'name' column should only contain strings in column data of experiment '" * key * "'"))
        end
        if !allunique(curnames)
            throw(ErrorException("column names in experiment '" * key *"' must be unique"))
        end

        for n in curnames
            push!(all_names, n)
        end
    end
    return all_names
end

function check_sampledata(sampledata)
    if size(sampledata)[2] < 0 || names(sampledata)[1] != "name"
        throw(ErrorException("'sampledata' should contain at least 1 column named 'name'"))
    end

    samp_names = sampledata[!,1]
    if !isa(samp_names, AbstractVector{<:AbstractString})
        throw(ErrorException("'name' column of 'sampledata' should be a string vector"))
    end

    if !allunique(samp_names)
        throw(ErrorException("'name' column of 'sampledata' should contain unique values"))
    end
end

function check_samplemap(samplemap)
    expected = ["sample", "experiment", "colname"]
    if size(samplemap)[2] != 3 || names(samplemap) != expected
        throw(ErrorException("'samplemap' should contain columns 'sample', 'experiment' and 'colname'"))
    end

    for field in expected
        if !isa(samplemap[!,field], AbstractVector{<:AbstractString})
            throw(ErrorException("'" * field * "' column of 'samplemap' should be a string vector"))
        end
    end
end

"""
The `MultiAssayExperiment` class is a Bioconductor container for multimodal studies.
This is basically a list of `SummarizedExperiment` objects, each of which represents a particular experimental modality.
A mapping table specifies the relationships between the columns of each `SummarizedExperiment` and a conceptual "sample",
assuming that each sample has data for zero, one or multiple modalities.
A sample can be defined as anything from a cell line culture to an individual patient, depending on the context.

The central idea is to use the sample mapping to easily filter the `MultiAssayExperiment` based on the samples of interest.
For example, a user can call [`multifilter`](@ref) to only keep the columns of each `SummarizedExperiment` that correspond to desired samples via the sample mapping.
This facilitates coordination across multiple modalities without needing to manually subset each experiment.
We also store sample-level annotations in a sample data `DataFrame`,
where they can be easily attached to the `coldata` of a `SummarizedExperiment` for further analyses.

This implementation makes a few changes from the original Bioconductor implementation.
We do not consider the `MultiAssayExperiment` to contain any "columns", as this was unnecessarily confusing.
The previous `colData` field has thus been renamed to `sampledata`, to reflect the fact that we are operating on samples.
We are also much more relaxed about harmonization between the experiments, sample mapping, and sample data -
or more specifically, we don't harmonize at all, allowing greater flexibility in storage and manipulation.
"""
mutable struct MultiAssayExperiment
    experiments::DataStructures.OrderedDict{String,SummarizedExperiments.SummarizedExperiment}
    sampledata::DataFrames.DataFrame
    samplemap::DataFrames.DataFrame
    metadata::Dict{String,Any}

    @doc """
        MultiAssayExperiment()

    Creates an empty `MultiAssayExperiment` object.

    # Examples
    ```jldoctest
    julia> using MultiAssayExperiments

    julia> MultiAssayExperiment()
    MultiAssayExperiment object
      experiments(0):
      sampledata(1): name
      metadata(0):
    ```
    """
    function MultiAssayExperiment()
        samp = DataFrames.DataFrame(sample = Vector{String}(), experiment = Vector{String}(), colname = Vector{String}())
        cd = DataFrames.DataFrame(name = Vector{String}())
        new(
            Dict{String,SummarizedExperiments.SummarizedExperiment}(),
            cd,
            samp,
            Dict{String,Any}()
        )
    end

    @doc """
        MultiAssayExperiment(experiments)

    Creates an `MultiAssayExperiment` object from a set of `experiments`.
    The per-sample column data and sample mapping is automatically created from the union of column names from all `experiments`.

    # Examples
    ```jldoctest
    julia> using MultiAssayExperiments

    julia> import SummarizedExperiments, DataStructures

    julia> exp = DataStructures.OrderedDict{String,SummarizedExperiments.SummarizedExperiment}();

    julia> exp["foo"] = SummarizedExperiments.exampleobject(100, 10);

    julia> exp["bar"] = SummarizedExperiments.exampleobject(50, 20);

    julia> out = MultiAssayExperiment(exp)
    MultiAssayExperiment object
      experiments(2): foo bar
      sampledata(1): name
      metadata(0):
    ```
    """
    function MultiAssayExperiment(experiments::DataStructures.OrderedDict{String,SummarizedExperiments.SummarizedExperiment})
        # Gather all the unique column names and use them to create a column data entry.
        all_names = harvest_all_colname(experiments)
        union = [x for x in all_names]
        dummy_cd = DataFrames.DataFrame(name = union)

        # Creating a dummy sample mapping.
        all_exp = Vector{String}()
        all_samp = Vector{String}()
        for (key, val) in experiments
            tmp = Vector{String}(undef, size(val)[2])
            fill!(tmp, key)
            append!(all_exp, tmp)

            cd = SummarizedExperiments.coldata(val)
            curnames = cd[!,"name"]
            append!(all_samp, curnames)
        end
        mapping = DataFrames.DataFrame(sample = all_samp, experiment = all_exp, colname = copy(all_samp))

        new(
            experiments,
            dummy_cd,
            mapping,
            Dict{String,Any}()
        )
    end

    @doc """
        MultiAssayExperiment(experiments, sampledata, samplemap, metadata = Dict{String,Any}())

    Creates a new `MultiAssayExperiment` from its components.

    `experiments` should contain ordered pairs of experiment names and `SummarizedExperiment` objects.
    Each `SummarizedExperiment` may contain any number and identity for the rows.
    However, the column names must be non-nothing and unique within each object.

    Each row of `sampledata` corresponds to a conceptual sample.
    The first column should be called `name` and contain the names of the samples in a `Vector{String}`.
    Sample names are arbitrary but should be unique.
    Any number and type of other columns may be provided, usually containing sample-level annotations.

    The `samplemap` table is expected to have 3 `Vector{String}` columns - `sample`, `experiment` and `colname` -
    specifying the correspondence between each conceptual sample and the columns of a particular `SummarizedExperiment`.
    See [`setsamplemap!`](@ref) for more details on the expected format.

    Note that values in the `samplemap` columns need not have a 1:1 match to their cross-referenced target;
    any values unique to one or the other will be ignored in methods like [`expandsampledata`](@ref) and [`filtersamplemap`](@ref).
    This allows users to flexibly manipulate the object without constantly hitting validity checks.

    The `metadata` stores other annotations unrelated to the samples.

    # Examples
    ```jldoctest
    julia> using MultiAssayExperiments

    julia> import SummarizedExperiments, DataStructures, DataFrames;

    julia> exp = DataStructures.OrderedDict{String,SummarizedExperiments.SummarizedExperiment}();

    julia> exp["foo"] = SummarizedExperiments.exampleobject(100, 2);

    julia> exp["bar"] = SummarizedExperiments.exampleobject(50, 5);

    julia> cd = DataFrames.DataFrame(
               name = ["Aaron", "Michael", "Jayaram", "Sebastien", "John"],
               disease = ["good", "bad", "good", "bad", "very bad"]
           );

    julia> sm = DataFrames.DataFrame(
               sample = ["Aaron", "Michael", "Aaron", "Michael", "Jayaram", "Sebastien", "John"],
               experiment = ["foo", "foo", "bar", "bar", "bar", "bar", "bar"],
               colname = ["Patient1", "Patient2", "Patient1", "Patient2", "Patient3", "Patient4", "Patient5"]
           );

    julia> using MultiAssayExperiments;

    julia> out = MultiAssayExperiment(exp, cd, sm)
    MultiAssayExperiment object
      experiments(2): foo bar
      sampledata(2): name disease
      metadata(0):
    ```
    """
    function MultiAssayExperiment(
            experiments::DataStructures.OrderedDict{String,SummarizedExperiments.SummarizedExperiment},
            sampledata::DataFrames.DataFrame,
            samplemap::DataFrames.DataFrame,
            metadata::Dict{String,Any} = Dict{String,Any}()
        )

        # Running through sanity checks.
        check_sampledata(sampledata)
        harvest_all_colname(experiments)
        check_samplemap(samplemap)

        new(
            experiments,
            sampledata,
            samplemap,
            metadata
        )
    end
end


text:
# This file is a part of JuliaFEM.
# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md

using JuliaFEM
using SparseArrays, Test

function test_auxiliary_plane_transforms()
    nodes = Vector{Float64}[
        [0.0, 0.0, 0.0],
        [1.0, 0.0, 0.0],
        [0.0, 1.0, 0.0]]
    e1 = Tri3([1, 2, 3])
    # local coordinate system N, T1, T2 in node
    R = [0.0 1.0 0.0
         0.0 0.0 1.0
         1.0 0.0 0.0]
    e1["geometry"] = Vector{Float64}[nodes[1], nodes[2], nodes[3]]
    e1["normal-tangential coordinates"] = Matrix{Float64}[R, R, R]
    time::Real = 0.0
    x0, Q = create_auxiliary_plane(e1, time)
    @info("x0 = $x0")
    @info("Q = $Q")
    @test isapprox(x0, [1.0/3.0, 1.0/3.0, 0.0])
    @test isapprox(Q, R)
    p1 = Float64[1.0/3.0+0.1, 1.0/3.0+0.1, 1.0]
    p2 = project_point_to_auxiliary_plane(p1, x0, Q)
    @info("point in auxiliary plane p2 = $p2")
    @test isapprox(p2, [0.1, 0.1])
    theta = project_point_from_plane_to_surface(p2, x0, Q, e1, time)
    @info("theta = $theta")
    @test isapprox(theta[1], 0.0)
    X = e1("geometry", theta[2:3], time)
    @info("projected point = $X")
    @test isapprox(X, Float64[1.0/3.0+0.1, 1.0/3.0+0.1, 0.0])
end

function test_get_edge_intersections()
    # first case, two triangles
    S = [ 0.0 0.0; 3.0  0.0; 0.0 3.0]'
    M = [-1.0 1.0; 2.0 -0.5; 1.0 1.5]'
    P, n = get_edge_intersections(S, M)
    P_expected = [
         1.00 1.75 0.00 0.00
         0.00 0.00 0.50 1.25]
    n_expected = [
        1 1 0
        0 0 0
        1 0 1]
    @test isapprox(P, P_expected)
    @test isapprox(n, n_expected)

    # slave 4 vertices non-convex, master triangle
    S = [ 0.0 0.0; 2.5  0.0; 1.0 1.0; 0.0 2.0]'
    M = [-1.0 1.0; 2.0 -0.5; 1.0 1.5]'
    P, n = get_edge_intersections(S, M)
    P_expected = [
        1.0 1.75 1.375 0.60 0.00 0.00
        0.0 0.00 0.750 1.40 0.50 1.25]
    n_expected = [
        1 1 0
        0 1 0
        0 0 1
        1 0 1]
    @test isapprox(P, P_expected)
    @test isapprox(n, n_expected)

    # slave 3 triangle, master 4 vertices
    S = [ 0.0 0.0; 3.0  0.0; 0.0 3.0]'
    M = [-1.0 1.0; 2.0 -0.5; 1.0 1.5; -1.0 2.0]'
    P, n = get_edge_intersections(S, M)
    P_expected = [
        1.00 1.75 0.00 0.00
        0.00 0.00 0.50 1.75]
    n_expected = [
        1 1 0 0
        0 0 0 0
        1 0 1 0]
    @test isapprox(P, P_expected)
    @test isapprox(n, n_expected)
end

function test_get_points_inside_triangle()
    S = [0.0 0.0; 3.0 0.0; 0.0 3.0]'
    pts = [-1.0 1.0; 2.0 -0.5; 1.0 1.5; 0.5 1.5]'
    P = get_points_inside_triangle(S, pts)
    @test isapprox(P, [1.0 1.5; 0.5 1.5]')
end

function test_is_point_inside_convex_polygon()
    X = Vector{Float64}[[0.0, 0.0], [1.0, 0.0], [1.0, 1.0], [0.0, 1.0]]
    @test is_point_inside_convex_polygon([0.5, 0.5], X) == true
    @test is_point_inside_convex_polygon([1.0, 0.5], X) == true
    @test is_point_inside_convex_polygon([1.1, 0.5], X) == false
    @test is_point_inside_convex_polygon([1.0, 1.0], X) == true
    @test is_point_inside_convex_polygon([0.0, 0.3], X) == true
    @test is_point_inside_convex_polygon([0.0, -0.000001], X) == false
end

function test_polygon_clipping_easy()
    S = [0 0; 3 0; 0 3]'
    M = [-1 1; 2 -1/2; 2 2]'
    P, n = clip_polygon(S, M)
    @test isapprox(P, [0.0 0.5; 1.0 0.0; 2.0 0.0; 2.0 1.0; 1.25 1.75; 0.0 4/3]')
    @test isapprox(n, [1 0 1; 1 1 0; 0 1 1])
end

function test_polygon_clipping_no_clip()
    # no clipping at all
    S = [-0.125   0.125  0.125  -0.125
         -0.125  -0.125  0.125   0.125]
    M = [-0.291667  -0.625     -0.625  -0.291667
         -0.208333  -0.208333   0.125   0.125   ]
    P, n = clip_polygon(S, M)
    # FIXME: check better.
    @test isa(P, Void)
    @test isa(n, Void)

end

function test_calculate_polygon_centerpoint()
    P = [
        0.0  1.0  2.0  2.0  1.25  0.0
        0.5  0.0  0.0  1.0  1.75  1.33333]
    C = calculate_polygon_centerpoint(P)
    @info("Polygon centerpoint: $C")
    @test isapprox(C, [1.0397440690338993, 0.8047003412233396])
end

function test_assemble_3d_problem_tri3()
    nodes = Vector{Float64}[
        [0.0, 0.0, 0.0],
        [1.0, 0.0, 0.0],
        [0.0, 1.0, 0.0],
        [0.0, 0.0, 0.1],
        [1.0, 0.0, 0.1],
        [0.0, 1.0, 0.1]]
    mel = Tri3([4, 5, 6])
    mel["geometry"] = Vector{Float64}[nodes[4], nodes[5], nodes[6]]
    sel = Tri3([1, 2, 3])
    sel["geometry"] = Vector{Float64}[nodes[1], nodes[2], nodes[3]]
#   Rv = [0.0 1.0 0.0
#         0.0 0.0 1.0
#         1.0 0.0 0.0]
#   sel["normal-tangential coordinates"] = Matrix{Float64}[Rv, Rv, Rv]
    calculate_normal_tangential_coordinates!(sel, 0.0)
    sel["master elements"] = Element[mel]
    prob = MortarProblem("temperature", 1)

    push!(prob, sel)
    stiffness_matrix = full(assemble(prob, 0.0).stiffness_matrix)
    @info("stiffness matrix for this problem:\n$stiffness_matrix")
    M = D = 1/24*[2 1 1; 1 2 1; 1 1 2]
    B = [D -M]  # slave dofs are first in this.
    @info("expected matrix for this problem:\n$B")
    @test isapprox(stiffness_matrix, B)

    # rotate and translate surface and check that we are still having same results
    Rx(t) = [
        1.0 0.0     0.0
        0.0 cos(t) -sin(t)
        0.0 sin(t)  cos(t)]
    Ry(t) = [
        cos(t)  0.0 sin(t)
        0.0     1.0 0.0
        -sin(t) 0.0 cos(t)
        ]
    Rz(t) = [
        cos(t) -sin(t) 0.0
        sin(t)  cos(t) 0.0
        0.0     0.0    1.0]
    T = [1.0, 1.0, 1.0]
    tx = pi/3.0
    ty = pi/4.0
    tz = pi/5.0
    for node in nodes
        node[:] = Rz(tz)*Ry(ty)*Rx(tx)*node + T
    end
    calculate_normal_tangential_coordinates!(sel, 0.0)
    stiffness_matrix = full(assemble(prob, 0.0).stiffness_matrix)
    @info("sel midpnt: ", sel("geometry", [1/3, 1/3], 0.0))
    @info("nt basis: ", sel("normal-tangential coordinates", [1/3, 1/3], 0.0))
    @test isapprox(stiffness_matrix, B)

end

function test_assemble_3d_problem_quad4()
    @info("assemble 3d problem in quad4-quad4")
    nodes = Vector{Float64}[
        [0.0, 0.0, 0.0],
        [1.0, 0.0, 0.0],
        [1.0, 1.0, 0.0],
        [0.0, 1.0, 0.0],
        [0.0, 0.0, 0.1],
        [2.0, 0.0, 0.1],
        [2.0, 2.0, 0.1],
        [0.0, 2.0, 0.1]]

    mel = Quad4([5, 6, 7, 8])
    mel["geometry"] = Vector{Float64}[nodes[5], nodes[6], nodes[7], nodes[8]]
    sel = Quad4([1, 2, 3, 4])
    sel["geometry"] = Vector{Float64}[nodes[1], nodes[2], nodes[3], nodes[4]]
    calculate_normal_tangential_coordinates!(sel, 0.0)
    sel["master elements"] = Element[mel]
    prob = MortarProblem("temperature", 1)

    push!(prob, sel)
    stiffness_matrix = full(assemble(prob, 0.0).stiffness_matrix)*144
    D = [16 8 4 8; 8  16 8 4;  4 8 16 8;  8 4 8 16]
    M = [25 5 1 5; 20 10 2 4; 16 8  4 8; 20 4 2 10]
    B = [D -M]  # slave dofs are first in this.
    @info("expected matrix for this problem:")
    dump(round(B, 3))

    @info("stiffness matrix for this problem:")
    dump(round(stiffness_matrix, 3))
    @test isapprox(stiffness_matrix, B)

end

function test_assemble_3d_problem_quad4_2()
    @info("assemble 3d problem in quad4-quad4")
    nodes = Vector{Float64}[
        [0.0, 0.0, 0.0],
        [1/4, 0.0, 0.0],
        [1/4, 1/4, 0.0],
        [0.0, 1/4, 0.0],
        [0.0, 0.0, 0.0],
        [1/3, 0.0, 0.0],
        [1/3, 1/3, 0.0],
        [0.0, 1/3, 0.0]]

    mel = Quad4([5, 6, 7, 8])
    mel["geometry"] = Vector{Float64}[nodes[5], nodes[6], nodes[7], nodes[8]]
    sel = Quad4([1, 2, 3, 4])
    sel["geometry"] = Vector{Float64}[nodes[1], nodes[2], nodes[3], nodes[4]]
    calculate_normal_tangential_coordinates!(sel, 0.0)
    sel["master elements"] = Element[mel]
    prob = MortarProblem("temperature", 1)

    push!(prob, sel)
    stiffness_matrix = full(assemble(prob, 0.0).stiffness_matrix)*589824
    D = [
        4096 2048 1024 2048
        2048 4096 2048 1024
        1024 2048 4096 2048
        2048 1024 2048 4096
        ]

    M = [
        5184 1728  576 1728
        3456 3456 1152 1152
        2304 2304 2304 2304
        3456 1152 1152 3456
        ]
    B = [D -M]  # slave dofs are first in this.
    @info("expected matrix for this problem:")
    dump(round(B, 3))

    @info("stiffness matrix for this problem:")
    dump(round(stiffness_matrix, 3))
    @test isapprox(stiffness_matrix, B)

end

function test_assemble_3d_problem_quad4_3()
    @info("assemble 3d problem in quad4-quad4")
    a = 1/4
    b = 1/3
    nodes = Vector{Float64}[
        [2*a,   a, 0],
        [3*a,   a, 0],
        [3*a, 2*a, 0],
        [2*a, 2*a, 0],
        [  b,   0, 0],
        [2*b,   0, 0],
        [2*b,   b, 0],
        [  b,   b, 0]]

    mel = Quad4([5, 6, 7, 8])
    mel["geometry"] = Vector{Float64}[nodes[5], nodes[6], nodes[7], nodes[8]]
    sel = Quad4([1, 2, 3, 4])
    sel["geometry"] = Vector{Float64}[nodes[1], nodes[2], nodes[3], nodes[4]]
    calculate_normal_tangential_coordinates!(sel, 0.0)
    sel["master elements"] = Element[mel]
    prob = MortarProblem("temperature", 1)

    push!(prob, sel)
    stiffness_matrix = full(assemble(prob, 0.0).stiffness_matrix)*186624*9
    D = [
        7904 3040 560 1456
        3040 2432 448  560
         560  448 128  160
        1456  560 160  416
        ]

    M = [
        504 1224 7956 3276
        144  720 4680  936
         18   90  990  198
         63  153 1683  693
        ]

    B = [D -M]  # slave dofs are first in this.
    @info("expected matrix for this problem:")
    dump(round(B, 3))

    @info("stiffness matrix for this problem:")
    dump(round(stiffness_matrix, 3))
    @test isapprox(stiffness_matrix, B)

end

function test_3d_problem()
    nodes = Vector{Float64}[
        [0.0, 0.0, 0.0],
        [1.0, 0.0, 0.0],
        [1.0, 1.0, 0.0],
        [0.0, 1.0, 0.0],
        [0.0, 0.0, 0.5],
        [1.0, 0.0, 0.5],
        [1.0, 1.0, 0.5],
        [0.0, 1.0, 0.5],

        [0.0, 0.0, 0.5],
        [1.0, 0.0, 0.5],
        [1.0, 1.0, 0.5],
        [0.0, 1.0, 0.5],
        [0.0, 0.0, 1.0],
        [1.0, 0.0, 1.0],
        [1.0, 1.0, 1.0],
        [0.0, 1.0, 1.0],
    ]
    el1 = Hex8([1, 2, 3, 4, 5, 6, 7, 8])
    el2 = Hex8([9, 10, 11, 12, 13, 14, 15, 16])
    sym121 = Quad4([1, 2, 3, 4])
    sym131 = Quad4([1, 2, 6, 5])
    sym132 = Quad4([9, 10, 14, 13])
    sym231 = Quad4([4, 1, 5, 8])
    sym232 = Quad4([12, 9, 13, 16])
    force = Quad4([14, 15, 16, 13])
    l2u = Quad4([5, 6, 7, 8])
    u2l = Quad4([9, 10, 11, 12])
    elements = Element[el1, el2, sym121, sym131, sym132, sym231, sym232, force, l2u, u2l]
    update!(elements, "geometry", nodes)
    el1["youngs modulus"] = el2["youngs modulus"] = 900.0
    el1["poissons ratio"] = el2["poissons ratio"] = 0.25
    sym121["displacement 3"] = 0.0
    sym131["displacement 2"] = sym132["displacement 2"] = 0.0
    sym231["displacement 1"] = sym232["displacement 1"] = 0.0
    force["displacement traction force 3"] = -100.0
    l2u["master elements"] = Element[u2l]
    calculate_normal_tangential_coordinates!(l2u, 0.0)
    fb = LinearElasticityProblem("two elastic blocks")
    push!(fb, el1, el2, force)
    bc = DirichletProblem("symmetry boundaries", "displacement", 3)
    push!(bc, sym121, sym131, sym132, sym231, sym232)
    tie = MortarProblem("tie contact between bodies", "displacement", 3)
    push!(tie, l2u)
    solver = DirectSolver("solution of elasticity problem")
    push!(solver, fb)
    push!(solver, bc)
    push!(solver, tie)
    solver.nonlinear_problem = false
    solver.method = :UMFPACK
    solver()
    X = el2("geometry", [1.0, 1.0, 1.0], 0.0)
    u = el2("displacement", [1.0, 1.0, 1.0], 0.0)
    @info("displacement at $X = $u")
    @test isapprox(u, 1/36*[1, 1, -4])
end

#=
@testset "plane quad4 projector tests" begin
    a = 1/2
    b = 1/3
    nodes = Dict{Int64, Vector{Float64}}(
        1 => [0.0, 0.0, 0.0],
        2 => [1/2, 0.0, 0.0],
        3 => [1.0, 0.0, 0.0],
        4 => [0.0, 1.0, 0.0],
        5 => [1/2, 1.0, 0.0],
        6 => [1.0, 1.0, 0.0],

        7 => [0.0, 0.0, 0.0],
        8 => [1/3, 0.0, 0.0],
        9 => [2/3, 0.0, 0.0],
       10 => [1.0, 0.0, 0.0],
       11 => [0.0, 1/2, 0.0],
       12 => [1/3, 1/2, 0.0],
       13 => [2/3, 1/2, 0.0],
       14 => [1.0, 1/2, 0.0],
       15 => [0.0, 1.0, 0.0],
       16 => [1/3, 1.0, 0.0],
       17 => [2/3, 1.0, 0.0],
       18 => [1.0, 1.0, 0.0],
    )
    sel1 = Quad4([1, 2, 5, 4])
    sel2 = Quad4([2, 3, 6, 5])
    mel1 = Quad4([7, 8, 12, 11])
    mel2 = Quad4([8, 9, 13, 12])
    mel3 = Quad4([9, 10, 14, 13])
    mel4 = Quad4([11, 12, 16, 15])
    mel5 = Quad4([12, 13, 17, 16])
    mel6 = Quad4([13, 14, 18, 17])
    update(Element[sel1, sel2, mel1, mel2, mel3, mel4, mel5, mel6], "geometry", nodes)
    calculate_normal_tangential_coordinates!(sel1, 0.0)
    calculate_normal_tangential_coordinates!(sel2, 0.0)
    prob = MortarProblem("temperature", 1)
    push!(prob, sel1)
    push!(prob, sel2)

    sel1["master elements"] = [mel1, mel2, mel4, mel5]
    sel2["master elements"] = [mel2, mel3, mel5, mel6]
    stiffness_matrix = full(assemble(prob, 0.0).stiffness_matrix)*2592*6
    @info("interface matrix:")
    dump(round(stiffness_matrix, 3))
    B = [
        864  432   0 432  216   0 -420 -375  -15    0 -504 -450  -18    0  -84  -75   -3    0
        432 1728 432 216  864 216 -120 -690 -690 -120 -144 -828 -828 -144  -24 -138 -138  -24
          0  432 864   0  216 432    0  -15 -375 -420    0  -18 -450 -504    0   -3  -75  -84
        432  216   0 864  432   0  -84  -75   -3    0 -504 -450  -18    0 -420 -375  -15    0
        216  864 216 432 1728 432  -24 -138 -138  -24 -144 -828 -828 -144 -120 -690 -690 -120
          0  216 432   0  432 864    0   -3  -75  -84    0  -18 -450 -504    0  -15 -375 -420
        ]
    @info("expected interface matrix:")
    dump(round(B, 3))
    @test isapprox(stiffness_matrix, B)
end
=#

#= TODO: Fix test.
@testset "plane quad4 projector master 3x3 slave 2x2" begin
    a = 1/2
    b = 1/3
    nodes = Dict{Int64, Vector{Float64}}(
        1 => [0*a, 0*a, 0.0],
        2 => [1*a, 0*a, 0.0],
        3 => [2*a, 0*a, 0.0],
        4 => [0*a, 1*a, 0.0],
        5 => [1*a, 1*a, 0.0],
        6 => [2*a, 1*a, 0.0],
        7 => [0*a, 2*a, 0.0],
        8 => [1*a, 2*a, 0.0],
        9 => [2*a, 2*a, 0.0],

       10 => [0*b, 0*b, 0.0],
       11 => [1*b, 0*b, 0.0],
       12 => [2*b, 0*b, 0.0],
       13 => [3*b, 0*b, 0.0],
       14 => [0*b, 1*b, 0.0],
       15 => [1*b, 1*b, 0.0],
       16 => [2*b, 1*b, 0.0],
       17 => [3*b, 1*b, 0.0],
       18 => [0*b, 2*b, 0.0],
       19 => [1*b, 2*b, 0.0],
       20 => [2*b, 2*b, 0.0],
       21 => [3*b, 2*b, 0.0],
       22 => [0*b, 3*b, 0.0],
       23 => [1*b, 3*b, 0.0],
       24 => [2*b, 3*b, 0.0],
       25 => [3*b, 3*b, 0.0],
    )
    sel1 = Quad4([1, 2, 5, 4])
    sel2 = Quad4([2, 3, 6, 5])
    sel3 = Quad4([4, 5, 8, 7])
    sel4 = Quad4([5, 6, 9, 8])
    mel1 = Quad4([10, 11, 15, 14])
    mel2 = Quad4([11, 12, 16, 15])
    mel3 = Quad4([12, 13, 17, 16])
    mel4 = Quad4([14, 15, 19, 18])
    mel5 = Quad4([15, 16, 20, 19])
    mel6 = Quad4([16, 17, 21, 20])
    mel7 = Quad4([18, 19, 23, 22])
    mel8 = Quad4([19, 20, 24, 23])
    mel9 = Quad4([20, 21, 25, 24])
    update!(Element[sel1, sel2, sel3, sel4, mel1, mel2, mel3,
                    mel4, mel5, mel6, mel7, mel8, mel9], "geometry", nodes)
    calculate_normal_tangential_coordinates!(sel1, 0.0)
    calculate_normal_tangential_coordinates!(sel2, 0.0)
    calculate_normal_tangential_coordinates!(sel3, 0.0)
    calculate_normal_tangential_coordinates!(sel4, 0.0)
    prob = MortarProblem("temperature", 1)
    push!(prob, sel1)
    push!(prob, sel2)
    push!(prob, sel3)
    push!(prob, sel4)

    master_elements = [mel1, mel2, mel3, mel4, mel5, mel6, mel7, mel8, mel9]
    sel1["master elements"] = master_elements
    sel2["master elements"] = master_elements
    sel3["master elements"] = master_elements
    sel4["master elements"] = master_elements
    B = sparse(assemble(prob, 0.0).stiffness_matrix, 25, 25)*46656
    B = full(B)
    D = B[1:9,1:9]
    M = B[1:9,10:end]
    @info("interface matrix D:")
    dump(round(D, 3))
    @info("interface matrix M:")
    dump(round(M, 3))
    D_expected = [
        1296  648    0  648  324    0    0    0    0
         648 2592  648  324 1296  324    0    0    0
           0  648 1296    0  324  648    0    0    0
         648  324    0 2592 1296    0  648  324    0
         324 1296  324 1296 5184 1296  324 1296  324
           0  324  648    0 1296 2592    0  324  648
           0    0    0  648  324    0 1296  648    0
           0    0    0  324 1296  324  648 2592  648
           0    0    0    0  324  648    0  648 1296]

    M_expected = [
        -784  -700   -28    0  -700  -625   -25     0   -28   -25    -1     0    0     0     0    0
        -224 -1288 -1288 -224  -200 -1150 -1150  -200    -8   -46   -46    -8    0     0     0    0
           0   -28  -700 -784     0   -25  -625  -700     0    -1   -25   -28    0     0     0    0
        -224  -200    -8    0 -1288 -1150   -46     0 -1288 -1150   -46     0 -224  -200    -8    0
         -64  -368  -368  -64  -368 -2116 -2116  -368  -368 -2116 -2116  -368  -64  -368  -368  -64
           0    -8  -200 -224     0   -46 -1150 -1288     0   -46 -1150 -1288    0    -8  -200 -224
           0     0     0    0   -28   -25    -1     0  -700  -625   -25     0 -784  -700   -28    0
           0     0     0    0    -8   -46   -46    -8  -200 -1150 -1150  -200 -224 -1288 -1288 -224
           0     0     0    0     0    -1   -25   -28     0   -25  -625  -700    0   -28  -700 -784]

    @info("D - D_expected")
    dump(D - D_expected)
    @info("M - M_expected")
    dump(M - M_expected)

    @test isapprox(D, D_expected)
    @test isapprox(M, M_expected)
    #=
    B = [
        864  432   0 432  216   0 -420 -375  -15    0 -504 -450  -18    0  -84  -75   -3    0
        432 1728 432 216  864 216 -120 -690 -690 -120 -144 -828 -828 -144  -24 -138 -138  -24
          0  432 864   0  216 432    0  -15 -375 -420    0  -18 -450 -504    0   -3  -75  -84
        432  216   0 864  432   0  -84  -75   -3    0 -504 -450  -18    0 -420 -375  -15    0
        216  864 216 432 1728 432  -24 -138 -138  -24 -144 -828 -828 -144 -120 -690 -690 -120
          0  216 432   0  432 864    0   -3  -75  -84    0  -18 -450 -504    0  -15 -375 -420
        ]
    @info("expected interface matrix:")
    dump(round(B, 3))
    @test isapprox(stiffness_matrix, B)
=#
end
=#


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

@doc raw"""
    ObjsBotProfileIcons(;
        image_36=nothing,
        image_48=nothing,
        image_72=nothing,
    )

    - image_36::String
    - image_48::String
    - image_72::String
"""
mutable struct ObjsBotProfileIcons <: SwaggerModel
    image_36::Any # spec type: Union{ Nothing, String } # spec name: image_36
    image_48::Any # spec type: Union{ Nothing, String } # spec name: image_48
    image_72::Any # spec type: Union{ Nothing, String } # spec name: image_72

    function ObjsBotProfileIcons(;image_36=nothing, image_48=nothing, image_72=nothing)
        o = new()
        validate_property(ObjsBotProfileIcons, Symbol("image_36"), image_36)
        setfield!(o, Symbol("image_36"), image_36)
        validate_property(ObjsBotProfileIcons, Symbol("image_48"), image_48)
        setfield!(o, Symbol("image_48"), image_48)
        validate_property(ObjsBotProfileIcons, Symbol("image_72"), image_72)
        setfield!(o, Symbol("image_72"), image_72)
        o
    end
end # type ObjsBotProfileIcons

const _property_map_ObjsBotProfileIcons = Dict{Symbol,Symbol}(Symbol("image_36")=>Symbol("image_36"), Symbol("image_48")=>Symbol("image_48"), Symbol("image_72")=>Symbol("image_72"))
const _property_types_ObjsBotProfileIcons = Dict{Symbol,String}(Symbol("image_36")=>"String", Symbol("image_48")=>"String", Symbol("image_72")=>"String")
Base.propertynames(::Type{ ObjsBotProfileIcons }) = collect(keys(_property_map_ObjsBotProfileIcons))
Swagger.property_type(::Type{ ObjsBotProfileIcons }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_ObjsBotProfileIcons[name]))}
Swagger.field_name(::Type{ ObjsBotProfileIcons }, property_name::Symbol) =  _property_map_ObjsBotProfileIcons[property_name]

function check_required(o::ObjsBotProfileIcons)
    (getproperty(o, Symbol("image_36")) === nothing) && (return false)
    (getproperty(o, Symbol("image_48")) === nothing) && (return false)
    (getproperty(o, Symbol("image_72")) === nothing) && (return false)
    true
end

function validate_property(::Type{ ObjsBotProfileIcons }, name::Symbol, val)
end