Building a metacommunity

Everything this package computes starts from a Metacommunity, and it is made of three things:

• types — what the individuals are, and how similar they are to each other;
• a partition — how the whole is divided into subcommunities;
• abundances — how much of each type is in each subcommunity.

If you supply only abundances, the other two are inferred: types become UniqueTypes (all wholly distinct from one another), and the partition becomes one subcommunity per column.

You have X, write Y

julia> using Diversity
julia> using LinearAlgebra
julia> counts = [10 0 0; 5 5 0; 0 5 10]          # 3 types down, 3 subcommunities across3×3 Matrix{Int64}:
 10  0   0
  5  5   0
  0  5  10
julia> Z = [1.0 0.5 0.0; 0.5 1.0 0.5; 0.0 0.5 1.0]3×3 Matrix{Float64}:
 1.0  0.5  0.0
 0.5  1.0  0.5
 0.0  0.5  1.0

julia> Metacommunity(counts)Metacommunity{Float64, Matrix{Int64}, Matrix{Float64}, UniqueTypes, Subcommunities} with 3 species in 3 subcommunities measuring Unique diversity.

Species names:
1, 2, 3

Subcommunity names:
1, 2, 3
julia> Metacommunity(counts, Z)Metacommunity{Float64, Matrix{Int64}, Matrix{Float64}, GeneralTypes{Float64, Matrix{Float64}, Vector{Int64}}, Subcommunities} with 3 species in 3 subcommunities measuring Arbitrary Z diversity.

Species names:
1, 2, 3

Subcommunity names:
1, 2, 3
julia> Metacommunity(counts, UniqueTypes(["ash", "oak", "elm"]),
                     Subcommunities(["north", "middle", "south"]))Metacommunity{Float64, Matrix{Int64}, Matrix{Float64}, UniqueTypes, Subcommunities} with 3 species in 3 subcommunities measuring Unique diversity.

Species names:
ash, oak, elm

Subcommunity names:
north, middle, south

Naming your types and subcommunities costs one line and pays for itself the moment you read a results table, because the names appear in the type_name and partition_name columns.

Counts or proportions

You can pass either, but they are treated differently, and the difference is deliberate.

julia> getabundance(Metacommunity(counts))          # counts: normalised silently3×3 Matrix{Float64}:
 0.285714  0.0       0.0
 0.142857  0.142857  0.0
 0.0       0.142857  0.285714

Integer counts are counts, so they are divided by their total without comment. Floating-point abundances are assumed to be relative already, so if they do not sum to one you get a warning:

julia> percolumn = counts ./ sum(counts, dims = 1)   # each column sums to 1 - a common mistake3×3 Matrix{Float64}:
 0.666667  0.0  0.0
 0.333333  0.5  0.0
 0.0       0.5  1.0
julia> sum(percolumn)                                # ... so the whole thing sums to 3, not 13.0

Handing that to Metacommunity gets you a warning and a correction:

julia> Metacommunity(percolumn);
┌ Warning: Abundances not normalised to 1, correcting...
└ @ Diversity Metacommunity.jl:124

Abundances here are relative to the whole metacommunity, not to each subcommunity — the column sums are the subcommunity weights, and normalising per column throws them away, making every subcommunity look the same size. If you see it, you should probably go back to the original data and normalise over the whole metacommunity.

A worked example

Ten sites and nine species, built so you can predict what each measure will say:

• Sites 1–8 lie along a gradient, each holding two or three species, turning over as you move along it.
• Reference holds all eight gradient species, evenly.
• Refuge is tiny and holds a single species found nowhere else.

julia> species = ["Sp $c" for c in "ABCDEFGH"];
julia> push!(species, "Refuge sp");
julia> sites = ["Site $i" for i in 1:8];
julia> push!(sites, "Reference"); push!(sites, "Refuge");
julia> abundance = zeros(Int, 9, 10);
julia> for j in 1:8, i in max(1, j - 1):min(8, j + 1)
           abundance[i, j] = 10                      # a moving window along the gradient
       end
julia> abundance[1:8, 9] .= 6;                       # Reference: every species, evenly
julia> abundance[9, 10] = 4;                         # Refuge: one species, found nowhere else
julia> mc = Metacommunity(abundance, UniqueTypes(species), Subcommunities(sites))Metacommunity{Float64, Matrix{Int64}, Matrix{Float64}, UniqueTypes, Subcommunities} with 9 species in 10 subcommunities measuring Unique diversity.

Species names:
Sp A, Sp B, Sp C...Sp H, Refuge sp

Subcommunity names:
Site 1, Site 2, Site 3...Reference, Refuge

That structure is much easier to see than to read off the loop:

using Plots
heatmap(sites, species, abundance, yflip = true, xrotation = 45, c = :Blues,
        title = "Individuals per species per site", colorbar_title = "count")

The dark band down the diagonal is the gradient: each site shares species with its neighbours and none at all with the far end. Reference is the one column touching every gradient species, at lower abundance. Refuge is the single cell in a row of its own — nothing else is in that site, and that species is nowhere else.

Before running anything, decide what you expect. Reference should be the most diverse in isolation and the most representative; Refuge should be the least diverse in isolation but the most distinctive, and should contribute most per individual. Now check:

julia> using DataFrames
julia> results = DataFrame(site = sites,
                           weight = round.(getweight(mc), digits = 3),
                           ᾱ = round.(norm_sub_alpha(mc, 1).diversity, digits = 2),
                           ρ̄ = round.(norm_sub_rho(mc, 1).diversity, digits = 2),
                           β = round.(raw_sub_beta(mc, 1).diversity, digits = 2),
                           γ = round.(sub_gamma(mc, 1).diversity, digits = 1))10×6 DataFrame
 Row │ site       weight   ᾱ        ρ̄        β        γ       
     │ String     Float64  Float64  Float64  Float64  Float64 
─────┼────────────────────────────────────────────────────────
   1 │ Site 1       0.074      2.0     0.22     0.33      8.9
   2 │ Site 2       0.11       3.0     0.36     0.31      8.4
   3 │ Site 3       0.11       3.0     0.4      0.28      7.6
   4 │ Site 4       0.11       3.0     0.4      0.28      7.6
   5 │ Site 5       0.11       3.0     0.4      0.28      7.6
   6 │ Site 6       0.11       3.0     0.4      0.28      7.6
   7 │ Site 7       0.11       3.0     0.36     0.31      8.4
   8 │ Site 8       0.074      2.0     0.22     0.33      8.9
   9 │ Reference    0.176      8.0     0.98     0.18      8.2
  10 │ Refuge       0.015      1.0     0.01     1.0      68.0

Two things worth reading off that table:

• Refuge scores β = 1, the maximum possible. Distinctiveness is 1 exactly when nothing outside the subcommunity resembles anything inside it — which is true here by construction. Its representativeness is correspondingly at its minimum, which is the subcommunity's own weight.
• Refuge has the lowest ᾱ and much the highest γ. A single-species site is as dull as a site can be in isolation, yet each of its individuals contributes far more to the diversity of the whole than an individual from anywhere else. Alpha and beta alone would have told you to ignore it.

How the partition changes the answer

The same abundances, divided three ways.

Undivided. With one subcommunity there is nothing to be distinct from, so every beta measure collapses to 1 and alpha equals gamma:

julia> undivided = Metacommunity(vec(sum(abundance, dims = 2)), UniqueTypes(species),
                                 Onecommunity())Metacommunity{Float64, Vector{Int64}, Matrix{Float64}, UniqueTypes, Onecommunity} with 9 species in 1 subcommunity measuring Unique diversity.

Species names:
Sp A, Sp B, Sp C...Sp H, Refuge sp

Subcommunity names:
1
julia> norm_sub_beta(undivided, 1).diversity1-element Vector{Float64}:
 1.0
julia> norm_sub_rho(undivided, 1).diversity1-element Vector{Float64}:
 1.0

Divided. The gamma diversity of the whole is unchanged — dividing a community does not alter what is in it — but the beta measures now carry information:

julia> meta_gamma(undivided, 1).diversity1-element Vector{Float64}:
 8.304589883293069
julia> meta_gamma(mc, 1).diversity1-element Vector{Float64}:
 8.304589883293062
julia> norm_meta_beta(mc, 1).diversity1-element Vector{Float64}:
 2.5115292759185666

Shattered. Now split Reference into two identical halves. This creates no new ecology: the two halves have the same composition as their parent, so an honest measure of "how many distinct subcommunities are there?" must not move.

julia> halves = hcat(abundance[:, 1:8], abundance[:, 9] .÷ 2, abundance[:, 9] .÷ 2,
                     abundance[:, 10])9×11 Matrix{Int64}:
 10  10   0   0   0   0   0   0  3  3  0
 10  10  10   0   0   0   0   0  3  3  0
  0  10  10  10   0   0   0   0  3  3  0
  0   0  10  10  10   0   0   0  3  3  0
  0   0   0  10  10  10   0   0  3  3  0
  0   0   0   0  10  10  10   0  3  3  0
  0   0   0   0   0  10  10  10  3  3  0
  0   0   0   0   0   0  10  10  3  3  0
  0   0   0   0   0   0   0   0  0  0  4
julia> shattered = Metacommunity(halves, UniqueTypes(species), Subcommunities(11))Metacommunity{Float64, Matrix{Int64}, Matrix{Float64}, UniqueTypes, Subcommunities} with 9 species in 11 subcommunities measuring Unique diversity.

Species names:
Sp A, Sp B, Sp C...Sp H, Refuge sp

Subcommunity names:
1, 2, 3...10, 11
julia> norm_meta_beta(mc, 1).diversity, norm_meta_beta(shattered, 1).diversity([2.5115292759185666], [2.5115292759185666])
julia> norm_meta_alpha(mc, 1).diversity, norm_meta_alpha(shattered, 1).diversity([3.3065869320817467], [3.3065869320817467])

The normalised measures do not move. The raw ones deliberately do:

julia> raw_meta_rho(mc, 1).diversity, raw_meta_rho(shattered, 1).diversity([3.6328925452090037], [4.1055892048434455])

This is the clearest way to see what "raw" and "normalised" mean. Raw redundancy should rise when you cut a subcommunity in two, because you have genuinely created two subcommunities that duplicate each other. Normalised measures answer the question "how many distinct subcommunities are there really?", and the answer is unchanged. The two differ by exactly the subcommunity's weight w: a raw measure keeps it and so sees size, a normalised one divides it out and cannot. Which you want depends on the question you are asking; see What this package does differently for why the framework guarantees the second.

Getting the results out

Every measure returns a long DataFrame in the same shape, whatever the measure and whatever the scale. That is convenient for combining results and inconvenient for reading, so to get back to a familiar types × subcommunities matrix, unstack it:

julia> ind = inddiv(NormalisedAlpha(mc), 1)90×8 DataFrame
 Row │ div_type  measure          q      type_level  type_name  partition_leve ⋯
     │ String    String           Int64  String      String     String         ⋯
─────┼──────────────────────────────────────────────────────────────────────────
   1 │ Unique    NormalisedAlpha      1  type        Sp A       subcommunity   ⋯
   2 │ Unique    NormalisedAlpha      1  type        Sp B       subcommunity
   3 │ Unique    NormalisedAlpha      1  type        Sp C       subcommunity
   4 │ Unique    NormalisedAlpha      1  type        Sp D       subcommunity
   5 │ Unique    NormalisedAlpha      1  type        Sp E       subcommunity   ⋯
   6 │ Unique    NormalisedAlpha      1  type        Sp F       subcommunity
   7 │ Unique    NormalisedAlpha      1  type        Sp G       subcommunity
   8 │ Unique    NormalisedAlpha      1  type        Sp H       subcommunity
  ⋮  │    ⋮             ⋮           ⋮        ⋮           ⋮             ⋮       ⋱
  84 │ Unique    NormalisedAlpha      1  type        Sp C       subcommunity   ⋯
  85 │ Unique    NormalisedAlpha      1  type        Sp D       subcommunity
  86 │ Unique    NormalisedAlpha      1  type        Sp E       subcommunity
  87 │ Unique    NormalisedAlpha      1  type        Sp F       subcommunity
  88 │ Unique    NormalisedAlpha      1  type        Sp G       subcommunity   ⋯
  89 │ Unique    NormalisedAlpha      1  type        Sp H       subcommunity
  90 │ Unique    NormalisedAlpha      1  type        Refuge sp  subcommunity
                                                   3 columns and 75 rows omitted
julia> unstack(ind, :type_name, :partition_name, :diversity)9×11 DataFrame
 Row │ type_name  Site 1    Site 2    Site 3    Site 4    Site 5    Site 6     ⋯
     │ String     Float64?  Float64?  Float64?  Float64?  Float64?  Float64?   ⋯
─────┼──────────────────────────────────────────────────────────────────────────
   1 │ Sp A            2.0       3.0     Inf       Inf       Inf       Inf     ⋯
   2 │ Sp B            2.0       3.0       3.0     Inf       Inf       Inf
   3 │ Sp C          Inf         3.0       3.0       3.0     Inf       Inf
   4 │ Sp D          Inf       Inf         3.0       3.0       3.0     Inf
   5 │ Sp E          Inf       Inf       Inf         3.0       3.0       3.0   ⋯
   6 │ Sp F          Inf       Inf       Inf       Inf         3.0       3.0
   7 │ Sp G          Inf       Inf       Inf       Inf       Inf         3.0
   8 │ Sp H          Inf       Inf       Inf       Inf       Inf       Inf
   9 │ Refuge sp     Inf       Inf       Inf       Inf       Inf       Inf     ⋯
                                                               4 columns omitted

The columns that identify a row are div_type, measure, q, type_level, type_name, partition_level and partition_name; diversity holds the answer. Selecting several orders at once and filtering afterwards is usually faster than repeated calls, because the individual diversities are computed once when the measure is built:

julia> profile = subdiv(NormalisedAlpha(mc), [0, 1, 2]);
julia> filter(row -> row.partition_name == "Reference", profile)3×8 DataFrame
 Row │ div_type  measure          q      type_level  type_name  partition_leve ⋯
     │ String    String           Int64  String      String     String         ⋯
─────┼──────────────────────────────────────────────────────────────────────────
   1 │ Unique    NormalisedAlpha      0  types                  subcommunity   ⋯
   2 │ Unique    NormalisedAlpha      1  types                  subcommunity
   3 │ Unique    NormalisedAlpha      2  types                  subcommunity
                                                               3 columns omitted

Now on real data

The dataset above was built so you could check the measures against your own expectations. Real data is the other way round — you do not know the answer, which is the point of measuring. Here is the same analysis on the European amphibian distributions that ship with SpatialEcology: 73 species across 1010 grid cells.

Any EcoBase assemblage can be measured directly, with no conversion:

julia> using Diversity, SpatialEcology, CSV, DataFrames
julia> file = joinpath(dirname(pathof(SpatialEcology)), "..", "data", "amph_Europe.csv");
julia> raw = CSV.read(file, DataFrame);
julia> amph = Assemblage(raw[!, 4:end], raw[!, 1:3], sitecolumns = false)[ Info: Matrix data assumed to be presence-absence
Assemblage with 73 species in 1010 sites

Species names:
Salamandra_salamandra, _Calotriton_asper, _Calotriton_arnoldi...Chioglossa_lusitanica, Pleurodeles_waltl

Site names:
1, 2, 3...1009, 1010
julia> countsubcommunities(amph), counttypes(amph)(1010, 73)

julia> rho = norm_sub_rho(amph, 1)[!, :diversity];
julia> gamma = sub_gamma(amph, 1)[!, :diversity];
julia> extrema(rho)(0.0015166819680988435, 0.7559785032793942)
julia> extrema(gamma)(10.648018648018521, 1369.918374086275)

Representativeness runs from 0.0015 to 0.76: some cells hold an assemblage much like Europe's as a whole, others almost nothing like it. Contribution per individual spans two orders of magnitude — the cells at the top are those holding species found almost nowhere else, which is exactly the "hidden diversity" the Refuge site showed in miniature.

Because the assemblage carries its own coordinates, the results can go straight onto a map:

using Plots
plot(norm_sub_rho(amph, 1), amph, title = "Representativeness of European amphibian faunas",
     markersize = 2)

Read that as a question about reserve selection: the dark cells are the ones whose amphibian fauna is least like Europe's overall. They are not necessarily the most species-rich — that is what norm_sub_alpha would show — and the difference between the two maps is precisely what the framework was built to expose.

Taking a view of part of one

view restricts a metacommunity to some of its types, some of its subcommunities, or both, and is how the rest of the EcoBase ecosystem subsets an assemblage. Indices, names and boolean masks all work:

julia> sub = view(mc, sites = ["Reference", "Refuge"])Diversity.SubAssemblage{Float64, UniqueTypes, Subcommunities, SubArray{Float64, 2, Matrix{Float64}, Tuple{Base.OneTo{Int64}, Vector{Int64}}, false}} with 9 species in 2 subcommunities measuring Unique diversity.

Species names:
Sp A, Sp B, Sp C...Sp H, Refuge sp

Subcommunity names:
Reference, Refuge
julia> norm_sub_alpha(sub, 1)[!, [:partition_name, :diversity]]2×2 DataFrame
 Row │ partition_name  diversity 
     │ String          Float64   
─────┼───────────────────────────
   1 │ Reference             8.0
   2 │ Refuge                1.0

The result is a genuine view: it aliases the original's abundances rather than copying them, and it holds only part of them, so they no longer sum to one. That is deliberate. Only a Metacommunity requires abundances summing to one, and the measures normalise what they read, so the subset is measured as a metacommunity in its own right. Convert it with Metacommunity(sub) if you want a cached object instead of a window.

• Diversity.API
• • Usage
  • The contract
• Diversity.Ecology
• • Usage
  • Faith's PD
• The framework
• • Origins
  • What this package does differently
  • What you supply
  • Ordinariness
  • The viewpoint parameter
  • The measures
  • What the measures show you
  • Properties
  • Special cases
• Genetic diversity
• • Usage
• Diversity.Hill
• • Usage
• Diversity.jl
• • A package for measuring and partitioning diversity
• Diversity.Jost
• • Usage
• Building a metacommunity
• • You have X, write Y
  • Counts or proportions
  • A worked example
  • How the partition changes the answer
  • Getting the results out
  • Now on real data
  • Taking a view of part of one
• Phylogenetic diversity
• • Usage
  • Raw and processed types
  • Faith's PD
• Coming from vegan
• • Before the table: three differences that matter more
  • Alpha diversity of single communities
  • Dissimilarity between two communities
  • Partitioning across many subcommunities
  • Similarity between types
  • Data preparation

• Diversity.API
• Diversity.Diversity
• Diversity.Ecology
• Diversity.Hill
• Diversity.Jost
• Diversity.ShortNames
• Diversity.individualDiversity
• Diversity.metacommunityDiversity
• Diversity.subcommunityDiversity
• Diversity.API.AbstractMetacommunity
• Diversity.API.AbstractPartition
• Diversity.API.AbstractTypes
• Diversity.AbstractGenetic
• Diversity.AbstractPhyloTypes
• Diversity.DiversityLevel
• Diversity.DiversityMeasure
• Diversity.Gamma
• Diversity.GeneralTypes
• Diversity.GeneralTypes
• Diversity.Metacommunity
• Diversity.NormalisedAlpha
• Diversity.NormalisedBeta
• Diversity.NormalisedRho
• Diversity.Onecommunity
• Diversity.PhyloBranches
• Diversity.PowerMeanMeasure
• Diversity.RawAlpha
• Diversity.RawBeta
• Diversity.RawRho
• Diversity.RelativeEntropyMeasure
• Diversity.Species
• Diversity.SubAssemblage
• Diversity.Subcommunities
• Diversity.Taxonomy
• Diversity.UniqueTypes
• Diversity.API._addedoutputcols
• Diversity.API._calcabundance
• Diversity.API._calcordinariness
• Diversity.API._calcsimilarity
• Diversity.API._countsubcommunities
• Diversity.API._counttypes
• Diversity.API._getabundance
• Diversity.API._getaddedoutput
• Diversity.API._getdiversityname
• Diversity.API._getmetaabundance
• Diversity.API._getmetaordinariness!
• Diversity.API._getordinariness!
• Diversity.API._getpartition
• Diversity.API._getscale
• Diversity.API._getsubcommunitynames
• Diversity.API._gettypenames
• Diversity.API._gettypes
• Diversity.API._getweight
• Diversity.API._subsetpartition
• Diversity.API._subsettypes
• Diversity.API.floattypes
• Diversity.API.mcmatch
• Diversity.API.typematch
• Diversity.Ecology.faith_pd
• Diversity.Ecology.generalisedfaith_pd
• Diversity.Ecology.generalisedjaccard
• Diversity.Ecology.generalisedrichness
• Diversity.Ecology.generalisedshannon
• Diversity.Ecology.generalisedsimpson
• Diversity.Ecology.gower
• Diversity.Ecology.jaccard
• Diversity.Ecology.pielou
• Diversity.Ecology.richness
• Diversity.Ecology.shannon
• Diversity.Ecology.simpson
• Diversity.GeneticType
• Diversity.Hill.hillnumber
• Diversity.Jost.jostalpha
• Diversity.Jost.jostbeta
• Diversity._getmeta
• Diversity.addedoutputcols
• Diversity.calcsimilarity
• Diversity.countsubcommunities
• Diversity.counttypes
• Diversity.diversity
• Diversity.getASCIIName
• Diversity.getFullName
• Diversity.getName
• Diversity.getabundance
• Diversity.getaddedoutput
• Diversity.getdiversityname
• Diversity.getmetaabundance
• Diversity.getmetaordinariness!
• Diversity.getordinariness!
• Diversity.getpartition
• Diversity.getsubcommunitynames
• Diversity.gettypenames
• Diversity.gettypes
• Diversity.getweight
• Diversity.hassimilarity
• Diversity.inddiv
• Diversity.meta_gamma
• Diversity.metadiv
• Diversity.norm_meta_alpha
• Diversity.norm_meta_beta
• Diversity.norm_meta_rho
• Diversity.norm_sub_alpha
• Diversity.norm_sub_beta
• Diversity.norm_sub_rho
• Diversity.powermean
• Diversity.qD
• Diversity.qDZ
• Diversity.raw_meta_alpha
• Diversity.raw_meta_beta
• Diversity.raw_meta_rho
• Diversity.raw_sub_alpha
• Diversity.raw_sub_beta
• Diversity.raw_sub_rho
• Diversity.sub_gamma
• Diversity.subdiv
• Diversity.vcf_dataframe