geoml.plots

The figures, in two backends that draw the same set under the same names: Explorer in matplotlib, to print, and Interactive in plotly, to look at. Both are built from a container plus a choice of variable, and several plotly figures can be linked on one page with Dashboard.

Three of them are only honest on data the model has not seen — accuracy, spread_check and variogram — since at a training location a model interpolates its own measurement.

Matplotlib

Looking at a data set before modelling it.

An Explorer is a choice of data and variables – one continuous (or vector) and one categorical – held on to, so that the figures can be asked for one after another without repeating it. The categorical variable is what splits and colours every other figure, which is the question worth asking of most geoscientific data: does this population behave as one, or as several?

eda = geoml.plots.Explorer(point, continuous="Elements",
                           categorical="Rock")
eda.histogram()
eda.pairs()
eda.pca(explained=0.9)
eda.scene()

These are matplotlib figures, meant to be saved and printed. Interactive is the same set of figures in plotly, for looking at on a screen.

class geoml.plots.explorer.Explorer(data, continuous=None, categorical=None, model=None, palette=None, cmap=None, uncertainty=None)

Bases: Selection

Exploratory figures for one data set and a choice of variables.

Takes its arguments from base.Selection: the container, a continuous and a categorical variable, a model for the figures that need one, and the colours to draw them in.

histogram(bins=25, figsize=None)

The distribution of the continuous variable, one panel per component.

Split by category when there is one: the populations are drawn over each other with a common set of bins, so their spreads can be compared rather than only their shapes.

Return type:

Figure

pairs(kind='scatter', alpha=0.7, bins=60, log_counts=False, log=False, principal_components=0, upper=None, figsize=None, size=6)

Every component against every other, with its distribution down the diagonal.

Parameters:

• kind (str) – “scatter”, or “hist2d” to count the points into cells instead, for when there are too many to draw one by one. Counts make a single surface, so the categories are pooled.

• alpha (float) – How transparent each point is, for kind=”scatter”.

• bins (int) – Cells along each axis, for kind=”hist2d”.

• log_counts (bool) – Colour the cells by the logarithm of the count. Worth turning on whenever a few cells hold most of the data, which is most of the time with a skewed variable: on a linear scale they take the whole colour range and everything else reads as empty.

• log (bool) – Draw the data on a log scale: the centred log-ratio for a composition, whose parts carry a constant sum, and an ordinary logarithm otherwise. Geochemical data is usually closer to symmetric this way, and it is what makes the principal components of a composition drawable – they are directions in the log-ratios, so this is the space they belong to.

• principal_components (int) – Draw this many principal components over the data, in the data’s own axes – the reverse of pca, which puts the data on the components’ axes. Each is a line through the mean, one standard deviation of that component long either way, so its length says how much of the spread it accounts for and its slant says which measurements move together. The sign of a component means nothing, hence a line rather than an arrow.

• upper (str) – What to put in the upper triangle, which is otherwise left empty: “hist2d” or “density” for where the mass is with the categories pooled, or “correlation” for the coefficient alone. The lower triangle says who each point is; the upper says where the data as a whole sits.

Return type:

Figure

pca(explained=0.9, log=False, kind='scatter', alpha=0.7, bins=60, log_counts=False, figsize=None, size=6)

The same pairs plot, on principal components instead of measurements.

Only the components carrying explained of the variance are drawn. Each panel also holds the loadings: an arrow per original column, showing what it contributes to the two components on the axes. A composition is opened up with the centred log-ratio first – see prepare.centred_log_ratio for why.

Parameters:

• explained (float) – Share of the total variance to reach, between 0 and 1.

• log (bool) – Take the components of the logarithms rather than of the measurements. A composition is opened up either way – there is no useful PCA of proportions, whose covariance the constant sum makes singular – so this decides the matter only for everything else.

• kind – As in pairs.

• alpha – As in pairs.

• bins – As in pairs.

• log_counts – As in pairs.

Return type:

Figure

scene(color=None, clip=None, figsize=None, size=14)

Where the data is, coloured by a variable.

The coordinates decide the drawing: a value against position in 1D, a map in 2D, a scatter in 3D. For anything more than a look at a 3D data set, as_pyvista() on the container is the better road.

Parameters:

• color (str) – The variable to colour by. Defaults to the continuous variable the Explorer holds, or the categorical one if that is all there is.

• clip (pair of floats) – Where to end the colour scale, as quantiles: [0, 0.99] for a variable with a long right tail, which is most assays. Without it one value far from the rest takes the whole scale and leaves everything else in a single shade. Nothing is dropped – the points beyond the ends take the end colour. Has no effect in 1D, where the value is an axis rather than a colour.

Return type:

Figure

training_curve(window=None, figsize=None)

The ELBO against the iteration it was measured at.

Every value is an estimate, from a sample of the latent variables and, under train_svi, a sample of the data too, so the curve is noisy whether or not training has settled. The running mean over it is the line to read: flat means finished, still climbing means it is not.

Parameters:

window (int) – Points to average over. Defaults to a fiftieth of the log.

Return type:

Figure

transformed_pairs(kind='scatter', alpha=0.7, bins=60, log_counts=False, upper=None, figsize=None, size=6)

The measurements as the model sees them, after its warping.

Two things worth checking before trusting a fitted model, and both are easier to see than to test. Down the diagonal, whether the warping made each variable Gaussian: the standard normal it is aiming at is drawn over each histogram. Off the diagonal, whether what is left is independent: a round cloud with a correlation near zero is what the model assumes, and a tilted or curved one is structure it will not capture.

The columns are numbered rather than named: a warping may rotate the data or bend it, so a column is generally a mixture of what was measured rather than any one of it.

Parameters:

• kind – As in pairs.

• alpha – As in pairs.

• bins – As in pairs.

• log_counts – As in pairs.

• upper (str) – What to put in the upper triangle, otherwise left empty: “hist2d”, “density” or “correlation”, as in pairs. “density” earns its place here: the contours of a pair the warping has done its work on are round and centred, and any lean or corner in them is the dependence the model is about to assume away.

Return type:

Figure

simulation_pairs(kind='hist2d', bins=60, log_counts=False, most=100000, margin=0.1, alpha=0.6, figsize=None, size=6)

What was measured against what was simulated.

The measurements fill the lower triangle and the simulations the upper, in the same form and between the same limits, so the two halves of the matrix can be read against each other: a simulation that reproduces the data has an upper half that mirrors the lower one. Down the diagonal the measured histogram carries the simulated density over it, which is the same comparison one variable at a time.

Both halves are counted into cells by default. Simulations come in location-by-realization blocks that run to millions of values, and a scatter of that many points is a filled rectangle whatever its transparency.

Parameters:

• kind (str) – “hist2d”, the default here, or “scatter”.

• bins – As in pairs.

• log_counts – As in pairs.

• alpha – As in pairs.

• size – As in pairs.

• most (int) – About how many simulated values to draw per component. They are taken by striding through the block, so the sample spans locations and realizations alike; the same stride is used for every component, since it is the pairs that are being looked at.

• margin (float) – How far past the measured range to look, as a share of it. Simulated values outside that window are left out: a few realizations reaching far beyond anything measured would otherwise set the scale for every panel and squeeze the comparison into a corner of it. They are dropped rather than pinned to the edge, which would pile the whole tail into the last cell and read as a mode that is not there.

Return type:

Figure

prediction_scatter(component=None, kind='scatter', alpha=0.6, bins=60, log_counts=False, figsize=None, size=10)

What was predicted against what was measured.

A single variable is drawn with its two distributions along the sides, which is where a bias shows that the scatter alone hides: the same cloud can sit on the 1:1 line while the predicted values are packed into a narrower range than the real ones. Several components are drawn as a panel each; name one with component to get the margins for it.

Parameters:

• component (str) – One component of a vector variable, drawn on its own.

• kind (str) – “scatter”, or “hist2d” for a block model, where the points are past counting and a scatter is a solid mass whatever its transparency.

• alpha (float) – How transparent each point is. Worth lowering as the points pile up, until there are enough to want kind=”hist2d” instead.

• bins (int) – Bins along each axis, for kind=”hist2d”.

• log_counts (bool) – Colour the cells by the logarithm of the count.

Return type:

Figure

accuracy(probabilities=None, figsize=None)

Whether the simulated spread is the spread the errors actually have.

For each probability, the share of true values that fall inside the interval holding that share of the simulations. On the 1:1 line the model knows what it does not know; below it the intervals are too narrow for the errors they have to cover, and above it the model is hedging. Deutsch’s goodness statistic sums that up in one number.

The intervals come from the model rather than from the container: what is stored there is the ground, the likelihood noise having been integrated out, and an assay is a measurement of the ground rather than the ground itself. Scoring the stored simulations against measured values would ask the model a question it never answered, and it would fail – so this figure needs the model the selection was built with.

Return type:

Figure

spread_check(bins=8, figsize=None)

Whether the noise the model fitted is the noise the data has.

A residual holds two things at once – how wrong the model was about the ground, and how far the assay fell from the ground – so it is read against the two together. The band is the noise, the line the whole claim, the points what the errors actually did. On the line is calibrated, below it is hedging, above it is over-confident.

The level axis is what says which term is at fault. A warping bends, so the noise grows with the value while the model’s own uncertainty does not: a shortfall widening with the grade is the noise, a flat one is the posterior. Points inside the band alone are the plainest case – the fitted noise over-explains the errors by itself.

Only honest on data the model has not seen: at a training location it interpolates its own measurement, and the residual is not an error.

Parameters:

bins (int or sequence) – How many bins, or where their edges are. A count gives equal-count bins, since a predicted grade is skewed and equal width would leave the top bins with a sample each; pass np.linspace(…) to ask for equal width instead.

Return type:

Figure

variogram(n_lags=15, max_lag=None, direction=None, tolerance=45.0, residuals=False, decluster=True, figsize=None)

The data’s spatial structure, against the fan the simulations make.

The experimental semivariogram of the measurements, with one thin curve per realization on the same pairs. A model that learned the spatial structure scatters its fan around the data’s curve; a kernel too smooth sags below it at short lags, and a nugget fitted into the range lifts it there. Neither shows in accuracy or spread_check, which judge one location at a time.

The measurements carry the likelihood noise and the realizations do not, so the fan is raised by what an independent error at each location adds to a semivariogram, taken from noise_variance. Without that the two curves are not the same quantity and every model looks over-smooth by a nugget.

With residuals=True it is the variogram of measured - predicted and the fan is dropped: structure left in the residuals is structure the model missed – honest on cross-validated predictions (models.cross_validate).

Each panel’s title carries VS, the variogram score of geoml.metrics.variogram_score() over the same locations and weights: the eye’s verdict on the fan as one number, for comparing two models without squinting. Lower is better, but only against another model on the same data – the score keeps a bias the curves are corrected for, and never reaches zero.

Parameters:

• n_lags (int) – Number of equal-width lag bins.

• max_lag (float, optional) – Longest separation considered; half the bounding-box diagonal by default.

• direction (array-like, optional) – Direction vector for a directional variogram; omnidirectional when absent. The anisotropy ellipsoid’s principal axes are the directions worth asking about.

• tolerance (float) – Angular tolerance around direction, in degrees.

• residuals (bool) – Variogram of the residuals instead, without the fan.

• decluster (bool or float) – Weight pairs by cell-declustering weights, so that the curve estimates the field’s variogram rather than the sampling’s. True chooses the cell size, a number fixes it, False leaves the pairs raw.

Return type:

Figure

grade_tonnage(component=None, density=None, cutoffs=30, max_uncertainty=None, log_mass=False, figsize=None)

How much material clears each cut-off, and how good it is.

Tonnage falls and grade rises as the cut-off climbs, and where the two cross is the question the curve is drawn to answer. Simulations are carried through one by one and drawn as a family, with the median over them picked out: the spread between the thin lines is what the model does not know about the answer.

Parameters:

• component (str) – Which grade, when the variable is a vector one. A cut-off applies to a single number, so there is nothing to guess here.

• density (float or str) – A number, a metadata column, or a ContinuousVariable – and in that last case its simulations are matched with the grade’s, one to one. Without a density the curve is in volume.

• cutoffs (int or array-like) – The grades to cut at, or how many of them to spread evenly across the range of the data.

• max_uncertainty (float) – Leave out the blocks the model doubts more than this, reading the column named when the Explorer was built. A block the model cannot speak for is not tonnage, and counting it flatters the answer at exactly the cut-offs where there is least data to go on.

• log_mass (bool) – Put the tonnage on a logarithmic scale. Most of a deposit clears the low cut-offs, so on a linear axis the high ones are a flat line along the bottom and the spread between the realizations there – which is where the decision usually is – cannot be seen at all. A cut-off that nothing clears has no logarithm and drops out of the curve rather than being drawn at the axis floor.

Return type:

Figure

Plotly

The same figures, for looking at rather than for printing.

Interactive answers every question Explorer does, with the same method names and the same arguments, and hands back a plotly figure instead of a matplotlib one:

eda = geoml.plots.Interactive(point, continuous="Elements",
                              categorical="Rock")
eda.pairs().show()
eda.scene()

What that buys is the things a printed figure cannot do: zooming a panel of a scatter matrix takes its whole row and column with it, so a cluster can be followed across every pair at once; clicking a category in the legend takes it out of every panel; and a three-dimensional scene can be turned around.

Every trace drawn from locations carries, in customdata, the row of the container each of its points came from. That is what makes a Dashboard able to link one figure’s selection to another’s – see dashboard. Figures whose points are not locations, such as the training curve, carry none, and a simulated value is a realization rather than a place, so the simulated half of simulation_pairs carries none either.

The numbers all come from prepare, the same functions Explorer reads, so the two backends cannot drift into showing different things.

class geoml.plots.interactive.Interactive(data, continuous=None, categorical=None, model=None, palette=None, cmap=None, uncertainty=None)

Bases: Selection

Exploratory figures for one data set, as plotly figures.

Takes its arguments from base.Selection: the container, a continuous and a categorical variable, a model for the figures that need one, and the colours to draw them in. Every method mirrors the one of the same name on Explorer, down to the arguments, except that a figure is sized in pixels (height, width) rather than in inches.

dashboard(figures=('scene', 'histogram', 'pairs'), **options)

The named figures on one page, sharing a selection.

Each name is a method of this class, called with its own defaults. That is the short way; for figures drawn with arguments of your own, or captioned, or from more than one data set, build a Dashboard directly – it takes figures rather than names:

geoml.plots.Dashboard(
    [("Where the cadmium is", eda.scene(color="Cd")),
     eda.pairs(log=True, upper="density")],
    title="Jura", columns=2)

Parameters:

• figures (sequence of str) – Which figures, in the order they are to be laid out. Asking for one this Interactive cannot draw raises whatever that figure would have raised on its own – there is no guessing here about what a data set can support.

• options – Passed to Dashboard: title, columns, plotlyjs, hint.

histogram(bins=25, height=None, width=None)

The distribution of the continuous variable, one panel per component.

Split by category when there is one: the populations are drawn over each other with a common set of bins, so their spreads can be compared rather than only their shapes.

Return type:

Figure

pairs(kind='scatter', alpha=0.7, bins=60, log_counts=False, log=False, principal_components=0, upper=None, size=5, height=None, width=None)

Every component against every other, with its distribution down the diagonal.

Parameters:

• kind (str) – “scatter”, or “hist2d” to count the points into cells instead, for when there are too many to draw one by one. Counts make a single surface, so the categories are pooled – and a counted cell is not a location, so a matrix drawn this way takes no part in a dashboard’s linked selection.

• alpha (float) – How opaque each point is, for kind=”scatter”.

• bins (int) – Cells along each axis, for kind=”hist2d”.

• log_counts (bool) – Colour the cells by the logarithm of the count. Worth turning on whenever a few cells hold most of the data, which is most of the time with a skewed variable.

• log (bool) – Draw the data on a log scale: the centred log-ratio for a composition, whose parts carry a constant sum, and an ordinary logarithm otherwise.

• principal_components (int) – Draw this many principal components over the data, in the data’s own axes – the reverse of pca. Each is a line through the mean, one standard deviation of that component long either way.

• upper (str) – What to put in the upper triangle, which is otherwise left empty: “hist2d” or “density” for where the mass is with the categories pooled, or “correlation” for the coefficient alone.

Return type:

Figure

pca(explained=0.9, log=False, kind='scatter', alpha=0.7, bins=60, log_counts=False, size=5, height=None, width=None)

The same pairs plot, on principal components instead of measurements.

Only the components carrying explained of the variance are drawn. Each panel also holds the loadings: an arrow per original column, showing what it contributes to the two components on the axes. A composition is opened up with the centred log-ratio first.

Parameters:

• explained (float) – Share of the total variance to reach, between 0 and 1.

• log (bool) – Take the components of the logarithms rather than of the measurements. A composition is opened up either way.

• kind – As in pairs.

• alpha – As in pairs.

• bins – As in pairs.

• log_counts – As in pairs.

Return type:

Figure

scene(color=None, clip=None, size=4, height=None, width=None)

Where the data is, coloured by a variable.

The coordinates decide the drawing: a value against position in 1D, a map in 2D, a scatter that can be turned around in 3D – which is the one plotly is worth reaching for on its own account.

A three-dimensional scene takes no part in a dashboard’s linked selection. Plotly has no box or lasso over a 3D scatter to select with, and no per-point selection state on one to show a selection made elsewhere. Several of them on one page do turn together, though.

Parameters:

• color (str or list) –

  The variable to colour by. Defaults to the continuous variable held, or the categorical one if that is all there is.

  A list of them puts a menu on the figure instead, one entry per choice, and the ground is drawn once with the values swapped over it. That is the answer to several variables over the same body of rock: five grades become one scene with five entries on a menu rather than five scenes, which is lighter, and it holds the point cloud still while the variable changes – the comparison one is actually making. A name in the list may be a variable, one component of a vector variable, or a vector variable standing for all of its components, so color=[“Elements”] names every grade in it.

  Only the locations measured in all of them are drawn; see prepare.color_choices.

• clip (pair of floats) – Where to end the colour scale, as quantiles: [0, 0.99] for a variable with a long right tail, which is most assays, and [0.01, 0.99] to take both ends. Without it one value far from the rest takes the whole scale and leaves everything else in a single shade. Nothing is dropped – the points beyond the ends take the end colour, and the hover still reports what was measured. A menu clips each of its choices by its own quantiles, since they are different variables and share nothing but the ground. Has no effect in 1D, where the value is an axis rather than a colour.

Return type:

Figure

training_curve(window=None, height=None, width=None)

The ELBO against the iteration it was measured at.

Every value is an estimate, from a sample of the latent variables and, under train_svi, a sample of the data too, so the curve is noisy whether or not training has settled. The running mean over it is the line to read: flat means finished, still climbing means it is not.

Parameters:

window (int) – Points to average over. Defaults to a fiftieth of the log.

Return type:

Figure

transformed_pairs(kind='scatter', alpha=0.7, bins=60, log_counts=False, upper=None, size=5, height=None, width=None)

The measurements as the model sees them, after its warping.

Two things worth checking before trusting a fitted model, and both are easier to see than to test. Down the diagonal, whether the warping made each variable Gaussian: the normal it is aiming at is drawn over each histogram. Off the diagonal, whether what is left is independent: a round cloud with a correlation near zero is what the model assumes, and a tilted or curved one is structure it will not capture.

The columns are numbered rather than named: a warping may rotate the data or bend it, so a column is generally a mixture of what was measured rather than any one of it.

Parameters:

• kind – As in pairs.

• alpha – As in pairs.

• bins – As in pairs.

• log_counts – As in pairs.

• upper – As in pairs.

Return type:

Figure

simulation_pairs(kind='hist2d', bins=60, log_counts=False, most=100000, margin=0.1, alpha=0.6, size=5, height=None, width=None)

What was measured against what was simulated.

The measurements fill the lower triangle and the simulations the upper, in the same form and between the same limits, so the two halves of the matrix can be read against each other: a simulation that reproduces the data has an upper half that mirrors the lower one. Down the diagonal the measured histogram carries the simulated density over it, which is the same comparison one variable at a time.

Both halves are counted into cells by default. Simulations come in location-by-realization blocks that run to millions of values, and a scatter of that many points is a filled rectangle whatever its transparency.

Parameters:

• kind (str) – “hist2d”, the default here, or “scatter”.

• bins – As in pairs.

• log_counts – As in pairs.

• alpha – As in pairs.

• size – As in pairs.

• most (int) – About how many simulated values to draw per component.

• margin (float) – How far past the measured range to look, as a share of it. Simulated values outside that window are left out.

Return type:

Figure

prediction_scatter(component=None, kind='scatter', alpha=0.6, bins=60, log_counts=False, size=6, height=None, width=None)

What was predicted against what was measured.

A single variable is drawn with its two distributions along the sides, which is where a bias shows that the scatter alone hides: the same cloud can sit on the 1:1 line while the predicted values are packed into a narrower range than the real ones. Several components are drawn as a panel each; name one with component to get the margins for it.

Parameters:

• component (str) – One component of a vector variable, drawn on its own.

• kind (str) – “scatter”, or “hist2d” for a block model, where the points are past counting.

• alpha (float) – How opaque each point is.

• bins (int) – Bins along each axis, for kind=”hist2d”.

• log_counts (bool) – Colour the cells by the logarithm of the count.

Return type:

Figure

accuracy(probabilities=None, height=None, width=None)

Whether the simulated spread is the spread the errors actually have.

For each probability, the share of true values that fall inside the interval holding that share of the simulations. On the 1:1 line the model knows what it does not know; below it the intervals are too narrow for the errors they have to cover, and above it the model is hedging. Deutsch’s goodness statistic sums that up in one number.

The intervals come from the model rather than from the container: what is stored there is the ground, the likelihood noise having been integrated out, and an assay is a measurement of the ground rather than the ground itself. So this figure needs the model the selection was built with.

Return type:

Figure

spread_check(bins=8, height=None, width=None)

Whether the noise the model fitted is the noise the data has.

A residual holds two things at once – how wrong the model was about the ground, and how far the assay fell from the ground – so it is read against the two together. The band is the noise, the line the whole claim, the points what the errors actually did. On the line is calibrated, below it is hedging, above it is over-confident.

The level axis is what says which term is at fault. A warping bends, so the noise grows with the value while the model’s own uncertainty does not: a shortfall widening with the grade is the noise, a flat one is the posterior.

Only honest on data the model has not seen.

Parameters:

bins (int or sequence) – How many bins, or where their edges are. A count gives equal-count bins; pass np.linspace(…) for equal width.

Return type:

Figure

variogram(n_lags=15, max_lag=None, direction=None, tolerance=45.0, residuals=False, decluster=True, height=None, width=None)

The data’s spatial structure, against the fan the simulations make.

The experimental semivariogram of the measurements, with one thin curve per realization on the same pairs. A model that learned the spatial structure scatters its fan around the data’s curve; a kernel too smooth sags below it at short lags, and a nugget fitted into the range lifts it there. Neither shows in accuracy or spread_check, which judge one location at a time.

The measurements carry the likelihood noise and the realizations do not, so the fan is raised by what an independent error at each location adds to a semivariogram, taken from noise_variance. Without that the two curves are not the same quantity and every model looks over-smooth by a nugget.

With residuals=True it is the variogram of measured - predicted and the fan is dropped: structure left in the residuals is structure the model missed – honest on cross-validated predictions (models.cross_validate).

Each panel’s title carries VS, the variogram score of geoml.metrics.variogram_score() over the same locations and weights: the eye’s verdict on the fan as one number, for comparing two models without squinting. Lower is better, but only against another model on the same data – the score keeps a bias the curves are corrected for, and never reaches zero.

Parameters:

• n_lags (int) – Number of equal-width lag bins.

• max_lag (float, optional) – Longest separation considered; half the bounding-box diagonal by default.

• direction (array-like, optional) – Direction vector for a directional variogram; omnidirectional when absent.

• tolerance (float) – Angular tolerance around direction, in degrees.

• residuals (bool) – Variogram of the residuals instead, without the fan.

• decluster (bool or float) – Weight pairs by cell-declustering weights, so that the curve estimates the field’s variogram rather than the sampling’s. True chooses the cell size, a number fixes it, False leaves the pairs raw.

Return type:

Figure

grade_tonnage(component=None, density=None, cutoffs=30, max_uncertainty=None, log_mass=False, height=None, width=None)

How much material clears each cut-off, and how good it is.

Tonnage falls and grade rises as the cut-off climbs, and where the two cross is the question the curve is drawn to answer. Simulations are carried through one by one and drawn as a family, with the median over them picked out: the spread between the thin lines is what the model does not know about the answer.

Two scales share one frame, so each gets gridlines in the colour of its own curve – a reader following the grade curve to a grey line would otherwise read a tonnage off it.

Parameters:

• component (str) – Which grade, when the variable is a vector one.

• density (float or str) – A number, a metadata column, or a ContinuousVariable – and in that last case its simulations are matched with the grade’s, one to one. Without a density the curve is in volume.

• cutoffs (int or array-like) – The grades to cut at, or how many of them to spread evenly across the range of the data.

• max_uncertainty (float) – Leave out the blocks the model doubts more than this, reading the column named when the Interactive was built.

• log_mass (bool) – Put the tonnage on a logarithmic scale. Most of a deposit clears the low cut-offs, so on a linear axis the high ones are a flat line along the bottom and the spread between the realizations there – which is where the decision usually is – cannot be seen. A cut-off that nothing clears has no logarithm and drops out.

Return type:

Figure

Several figures on one page, sharing a selection.

A figure answers one question. Most questions worth asking of spatial data are about two things at once – those samples, the ones out on the tail of the histogram, where are they on the map? – and no single figure answers that. A dashboard does: drag a box over any panel and the same locations light up in every other one.

eda = geoml.plots.Interactive(point, continuous="Elements",
                              categorical="Rock")
board = eda.dashboard()
board.write_html("jura.html")

board renders itself in a notebook, and write_html writes a page that carries plotly with it, so it can be sent to someone who has neither geoML nor a network.

Take the long way for control over what goes on it, and over how it is laid out – a nested entry is a row of its own, so the rows need not be equal:

geoml.plots.Dashboard(
    [("Where the samples are", eda.scene(color="Cd")),
     [eda.histogram(), eda.training_curve()],
     eda.pairs(log=True, upper="density")],
    title="Jura, before modelling")

What links the panels is customdata: every trace Interactive draws from locations carries the row of the container each of its points came from, and the script this module writes matches selections by those rows. Anything without them – a counted cell, a simulated value, a training curve, a 3D scene – is simply left alone when a selection is made, which is the honest answer for a panel whose points are not places.

3D scenes are linked the other way. There is no brushing them, but there is turning them, and a page of several answers a question one of them cannot: the same body of ground under two variables, seen from the same angle. Turn one and the rest follow.

class geoml.plots.dashboard.Dashboard(figures, title=None, columns=2, plotlyjs='embed', hint=True)

Bases: object

A page of figures whose selections are shared.

Parameters:

• figures (list) –

  The figures, in the order they are to be laid out. An entry may be:

  • a figure;

  • a `(caption, figure)` pair – the caption replaces whatever title the figure came with, which is how the same figure serves in one place under its own name and in another under the question it is there to answer;

  • a list or tuple of either, which is a row of its own, however many panels are in it.

  Nest one entry and the whole list is read as rows, a bare figure counting as a row of one; nest none and the figures are dealt out columns to a row. Rows are what a page usually wants: a map beside a histogram, and the scatter matrix that goes with them across the full width underneath.

  Dashboard([eda.scene(),

  [eda.histogram(), eda.training_curve()], eda.pairs()])

• title (str) – The heading over the page.

• columns (int) – How many panels to a row, when the rows were not given outright. A narrow window falls back to one, whatever this says.

• plotlyjs (str) – “embed” writes plotly into the page, some four megabytes of it, and the result opens anywhere with no network at all – which is what makes it a report rather than a link. “cdn” fetches the library instead, leaving a file small enough to keep in a repository or a notebook, and needing a connection to open.

• hint (str or bool) – The line under the heading saying what can be done with the page. False leaves it out.

height(page=1000)

How tall the page comes out, on one page pixels wide.

Only an estimate, and only wanted for the height to give the frame a notebook shows this in: the panels are laid out by the browser, and what they end up as depends on the width they are given. A row is as tall as the tallest panel in it, and a panel is as tall as its own proportions make it once it has been given its share of the width.

to_html(full=True)

The page, as text.

Parameters:

full (bool) – A whole document, head and all. False gives the piece that goes inside one, for pasting into a page of your own – and then plotly has to be loaded by that page, since a fragment carries no head to put it in.

write_html(path)

Write the page to a file, and hand back the path.

show(path=None)

Write the page and open it in a browser.

The arithmetic behind them

prepare holds what both backends read and imports no plotting library, so the numbers behind a figure can be had without drawing it — which is also how they are tested.

The numbers behind the figures.

Nothing here draws anything, and nothing here imports matplotlib. A figure is one way of looking at these arrays; a dashboard is another, and it reads the same functions rather than working the values out a second time. It also means the arithmetic – which components a fraction of variance asks for, what a composition looks like once it is opened up – can be tested against numbers instead of against pictures.

geoml.plots.prepare.variable(container, name)

The variable called name, or a message saying what there is.

geoml.plots.prepare.variable_or_component(container, name)

The variable called name, or the component of a vector variable.

A composition is held as one variable, so its parts are not in container.variables and asking for “Zn” would otherwise mean reaching into Elements by hand. The search is the container’s own – this module kept a near-copy of it for years, which is the duplication the path work was started to remove.

geoml.plots.prepare.numeric_values(var)

A continuous or vector variable as a matrix, and the rows that hold one.

Returns:

• values (array) – (n_data, n_columns), whatever the variable’s length.

• measured (array) – Boolean, (n_data,). A vector variable is measured at a location only where every component is.

• labels (list) – A name per column.

geoml.plots.prepare.category_values(var)

A categorical variable as one label per location.

Returns:

• values (array) – The label at each location; the empty string where nothing was measured, which is how a missing code reads.

• measured (array) – Boolean, (n_data,).

• labels (list) – The categories actually present, in the variable’s own order.

geoml.plots.prepare.groups(values, measured, labels)

A boolean mask per category, in the order the labels come in.

geoml.plots.prepare.centred_log_ratio(values)

A composition opened up into ordinary numbers.

Proportions carry a constant sum, so their covariance is negative by construction and a PCA of them describes the constraint as much as the data. The log-ratios have no such constraint. This is warping’s CenteredLogRatio in NumPy, for exploring rather than for modelling.

Return type:

NDArray[float64]

geoml.plots.prepare.logarithm(values, compositional=False)

The data on a log scale, by the road that suits it.

A composition takes the centred log-ratio: its parts carry a constant sum, and logging them one at a time leaves that constraint in place, whereas dividing each by the row’s geometric mean first removes it. Anything else takes an ordinary logarithm.

Non-positive values stop this rather than quietly becoming infinities: which columns they are in is worth knowing, since a zero in an assay usually means below detection rather than absent, and what to put in its place is a decision about the data, not about the figure.

Return type:

NDArray[float64]

geoml.plots.prepare.principal_components(values, explained=0.9)

The principal components of values, down to a share of the variance.

Follows warping.PCA: an eigendecomposition of the covariance, taken largest first. The scores are the plain projection, which is what a biplot puts on its axes; the loadings are the eigenvectors, for the arrows.

Parameters:

• values (array) – (n_data, n_columns).

• explained (float) – The share of the total variance to reach. The number of components is the fewest that reach it – never fewer than two, so that there is a plot to draw, and never more than there are columns.

Returns:

dict – scores (n_data, n_components), loadings (n_columns, n_components), ratio (the share each component carries, all of them), n_components, and – for drawing the components back onto the data they came from – mean, the centre they turn about, and eigenvalues, the variance along each, whose square root is a length in the data’s own units.

geoml.plots.prepare.component_analysis(var, explained=0.9, log=False)

principal_components for a variable, on a log scale if asked.

A composition is opened up whether or not log is set, and that is not an oversight: its parts carry a constant sum, so their covariance is singular and the first component of the raw proportions would spend itself describing the closure rather than the data. There is no useful PCA of a composition to refuse it. Everything else is left as measured unless log says otherwise.

geoml.plots.prepare.warped_values(model, name)

The measurements as the model sees them, after its warping.

A model does not work on the data as measured: the likelihood’s warping takes it somewhere it can be treated as Gaussian and, for several variables at once, uncorrelated. Whether it succeeded is a thing to look at rather than assume, and looking at it means passing the data through the same warping the model trained with – already fitted, since VGPNetwork initializes each likelihood from the measurements.

Only measured rows go through. get_measurements fills the rest with 1.0 to keep the array rectangular, and a warping has no reason to be kind to a value that was never there.

Returns:

• values (array) – (n_measured, warping.size_out).

• measured (array) – Boolean, (n_data,), saying which rows those are.

• labels (list) – A name per column.

geoml.plots.prepare.simulation_sample(var, most=100000)

The simulated values, thinned to a number that can be drawn.

A simulated variable holds n_data x n_sim values per component, which on a block model runs to hundreds of millions – more than a figure can show and, spread across every pair of a matrix, more than memory should be asked to hold. What a distribution looks like is settled by far fewer, so a stride is taken through the block – and taken while reading, a band of locations at a time, so the values that are thrown away are never all in memory together either.

The stride is the same for every component, and the mask that drops non-finite values is applied to whole rows. Thinning the components separately would pair the value simulated at one location with the value simulated at another, which is not a pair the model ever produced – the joint shape, which is the thing being looked at, would be an artefact.

Returns:

values (array) – (n_sample, n_columns), columns in the variable’s own order.

Return type:

NDArray[float64]

geoml.plots.prepare.padded_range(values, margin=0.1)

The span of each column, opened up by margin at both ends.

Returns a (low, high) pair per column. A column of one repeated value has no span to open up, so it is given a unit of room rather than a window of zero width.

geoml.plots.prepare.color_limits(values, clip=None)

The ends of a colour scale, set by where the data mostly is.

One value far from the rest takes the whole of a colour scale with it and leaves everything else in a single shade – the usual fate of a geochemical assay, whose tail is long and whose interest is not in it. Naming a pair of quantiles ends the scale where the data mostly ends instead: [0, 0.99] for a variable skewed to the right, [0.01, 0.99] to take both tails.

Nothing is dropped and nothing is altered. The values keep their own numbers and a hover still reports what was measured; what is bounded is the scale, so the few beyond it take the end colour rather than setting it.

Returns (low, high), or None when clip is None – and None again when the two ends come out equal, since a scale of zero width is no scale.

geoml.plots.prepare.color_choices(container, names)

Several variables’ values over the locations that carry all of them.

For a figure that keeps one set of points and swaps the values over them. The locations are those measured in every one of the variables named, which is the point: a cloud that gained and lost points as the choice changed would be two things changing at once, and the comparison – this variable against that one, here – is exactly what would be lost.

A name may be a variable, one component of a vector variable, or a vector variable itself, which stands for all of its components in order. So [“Elements”] names seven grades and [“Cd”, “Zn”] names two.

Returns:

• values (array) – (n_kept, n_columns), a column per choice in the order asked for.

• rows (array) – Which locations those are, as indices into the container.

• labels (list) – A name per column.

Return type:

tuple[NDArray[float64], NDArray[integer], list]

geoml.plots.prepare.cells(n_points, most)

How many cells to count n_points into, at most most of them.

The two halves of a comparison rarely hold comparable numbers – a few hundred measurements against a hundred thousand simulated values – and binning both the same way leaves the sparse half as scattered single counts that read as noise. Roughly the square root of the count keeps several points in a typical cell either way.

geoml.plots.prepare.counts_2d(x, y, bins=60, log=False)

Points counted into cells, with the empty ones left empty.

A cell holding nothing comes back as NaN rather than as zero, so that whatever draws it can leave it unpainted: painted, it takes the bottom of the colour scale and fills the panel with a background that looks like data.

Returns:

dict – x and y, the cell centres; z, what the colour is to be taken from, which is the count or its base-ten logarithm; and count, the count itself, so that a hover can say how many points a cell holds whichever of the two is being coloured. z and count are shaped (n_y, n_x), the way an image is indexed.

geoml.plots.prepare.density_grid(x, y, grid=60, most=4000)

Smoothed density over a pair of columns, on a regular mesh.

The estimate costs one pass over the data for every mesh cell, so a long column is thinned first: a density is a shape, and a few thousand points settle it as well as a few million do.

Returns (x_axis, y_axis, density) with density shaped (grid, grid) – rows along y_axis, as an image is – or None when there is nothing to smooth: fewer than three points, or a column that never varies, which has no width for a kernel to sit in.

geoml.plots.prepare.density_curve(values, limits, points=200, most=5000)

A smoothed distribution as a line, across limits.

Thinned as density_grid is, and for the same reason. Returns (grid, density), or None when there is nothing to smooth.

geoml.plots.prepare.normal_curve(values, low, high, points=200)

The normal of the same mean and spread as values, across a window.

Fitted rather than standard, so that the curve asks about the shape alone. A warping’s parameters are trained along with everything else and need not leave the data at unit variance – the GP’s amplitude absorbs a scale factor – so a standard normal would call a perfectly symmetric result skewed. Returns (x, density), or None for a column that never varies.

geoml.plots.prepare.continuous_parts(var)

The continuous columns a figure draws, one per component.

A scalar variable is its own single part; a vector or compositional one comes back as its components, in label order. Everything that reads measurements, prediction, noise_variance or the simulations asks through this, so a figure handed a categorical variable is refused by name rather than failing later at whichever column it reached for first.

Parameters:

var – A variable, as variable() returns it.

Returns:

list of geoml.data.ContinuousVariable – One part for a scalar variable, one per component otherwise.

Raises:

ValueError – If the variable carries no continuous columns of its own.

Return type:

list[ContinuousVariable]

geoml.plots.prepare.component_names(var)

The components a caller can name, ready to drop into a message.

Every figure that needs one grade out of several says so the same way, and each of them used to reach for labels – an attribute only some variables have. Asking through continuous_parts keeps that reach in one place.

Return type:

str

geoml.plots.prepare.prediction_values(container, name)

What was measured against what was predicted, component by component.

Returns:

• measured_values, predicted_values (array) – (n_compared, n_columns), holding only the locations that carry both.

• labels (list) – A name per column.

• rows (array) – Which locations those are, as indices into the container. A figure that can be brushed needs to say which row each of its points came from, and only the rows carrying both a measurement and a prediction are drawn.

geoml.plots.prepare.step_path(lo, hi, values)

A per-bin value as a polyline that steps at the edges.

One number per bin is not a curve, and drawing it through the bin centres says it is. Stepping at the edges shows where the bins are without a second thing on the figure to say so. A gap between bins comes back as a break rather than a line across it.

geoml.plots.prepare.spread_check(container, name, bins=8)

What a model claims a value’s spread is, against what it turned out to be.

A residual holds two things at once – how wrong the model was about the ground, and how far the assay fell from the ground – so it can only be read against the two together. This lays all three out along the predicted value: the noise the model fitted, the whole spread it claims, and the spread the errors actually had.

Reading it: the observed points on the claimed line means calibrated, below it means hedging, above it means over-confident. The level axis is what says which term is at fault. A warping bends, so the noise grows with the value while the model’s own uncertainty does not, and a shortfall that widens with the grade is the noise where a flat one is the posterior. Observed points sitting inside the noise band alone are the plainest case of all: the fitted noise over-explains the errors by itself.

Only honest on data the model has not seen. At a training location the model interpolates its own measurement and the residual is not an error.

Parameters:

• container (_SpatialData) – Point data carrying measurements and a prediction.

• name (str) – The variable.

• bins (int or sequence) – How many bins, or where their edges are. A count gives equal-count bins; positions are taken as they come, so np.linspace(…) is how to ask for equal width.

Returns:

list of dict – One per component, with label, the bin bounds lo/hi, the mean predicted value in each centre, the count, the observed root mean square residual and its observed_error, and the claimed noise and total spreads.

Return type:

list[dict]

geoml.plots.prepare.variogram(container, name, n_lags=15, max_lag=None, direction=None, tolerance=45.0, max_pairs=2000000, residuals=False, decluster=True)

The data’s spatial structure, against the fan the simulations reproduce.

The experimental semivariogram of the measured values, and one curve per realization computed on the same pairs: a model that learned the spatial structure scatters its realizations around the data’s curve, a kernel too smooth sags below it at short lags, and a nugget fitted into the range lifts it there. Neither shows in the marginal checks (accuracy, spread_check), which is what this figure exists for.

The two curves are put on the same footing before being compared. The measurements carry the likelihood noise and the stored realizations do not (predict integrates it out), so the fan is raised by what that noise adds to a semivariogram. For noise independent between locations that is exactly the pair-averaged (var_i + var_j) / 2, taken from the noise_variance column and added bin by bin – no draw, no seed, and exact in expectation rather than approximated. Without the correction the fan sits a nugget below the data at every lag and every model looks over-smooth. A container predicted without include_noise carries no such column and its fan is left where it is.

Pairs are declustered by default, which is the other half of making the comparison fair and matters more than it sounds. Samples follow the ore, so an experimental variogram computed on raw pairs describes the sampling as much as the field: on the bundled Walker Lake set, whose exhaustive truth is known, the raw sample curve runs 1.4 times the true variogram at the sill and 2.3 times at the shortest lag, and a model matching the field perfectly would look far too smooth against it. Each pair is therefore weighted by w_i * w_j from geoml.math.geometry.declustering_weights(), and the sill likewise, with the same weights used for the fan so that both sides estimate the same thing. Pass decluster=False for the raw curve, or a number to fix the cell size rather than let it be chosen.

With residuals=True the variogram is of measured - predicted instead and the fan is omitted: structure left in the residuals is structure the model missed, and on cross-validated predictions (models.cross_validate) it is honest. Whatever the flag, at a training location a model interpolates its own measurement, so the fan is only worth reading against data the model has not seen or with the data curve as the anchor.

Parameters:

• container (_SpatialData) – Point data carrying measurements (and simulations, for the fan).

• name (str) – The variable.

• n_lags (int) – Number of equal-width lag bins between zero and max_lag.

• max_lag (float, optional) – The longest separation considered. Half the bounding-box diagonal by default.

• direction (array-like, optional) – A direction vector for a directional variogram. Omnidirectional when absent. Call once per direction to compare them – the anisotropy ellipsoid’s principal axes are the ones worth asking about.

• tolerance (float) – Angular tolerance around direction, in degrees.

• max_pairs (int) – The pair budget. Past it the locations are strided down – deterministic, so two calls agree.

• residuals (bool) – Variogram of measured - predicted rather than of the measurements, with no fan.

• decluster (bool or float) – Weight pairs by cell-declustering weights, so that the curve estimates the field’s variogram rather than the sampling’s. True chooses the cell size, a number fixes it, False leaves the pairs raw.

Returns:

list of dict – One per component, with label, the bin centres lag, the pair count per bin, the data curve data, the sample variance sill, and realizations – an (n_realizations, n_lags) array, or None when there is nothing simulated or residuals was asked. noise is what was added to the fan per bin to put it on the measurements’ footing, or None when the container could not say. cell is the declustering cell used, or None when the pairs were left raw. score is geoml.metrics.variogram_score() over the same locations and weights, or None when there was no fan to score.

Return type:

list[dict]

Notes

The score is the figure’s verdict as one number, and it is a ranking rather than a measurement: unlike the curves it cannot be put on the measurements’ footing, since |difference| ** p is not a second moment and has no constant to add. It never reaches zero, and only comparisons between models on the same data mean anything. It is also taken over every pair of the locations kept, not only the binned ones, so it does not answer to max_lag or move when direction does.

geoml.plots.prepare.moving_average(values, window)

The running mean of values, and where each point belongs.

Returns the positions as well: a mean over window points only exists once there are that many, so the curve starts later than the one it smooths and has to be drawn against its own x.

Return type:

tuple[NDArray[integer], NDArray[float64]]

geoml.plots.prepare.training_curve(model, window=None)

The training log, with a running mean over it.

Each entry is one optimizer step – and under train_svi one batch, so the curve is noisy by construction: every value is the ELBO estimated from a sample of the data and a sample of the latent variables. The mean is what says whether it is still climbing.

Parameters:

window (int) – Points to average over. Defaults to a fiftieth of the log, which keeps the smoothing proportionate to however long training ran, but never fewer than five: a mean of two or three smooths nothing, and draws a second line on top of the first saying the same thing. A window longer than the log leaves it alone, and then there is no second line at all.

geoml.plots.prepare.realizations(var)

A variable’s simulations, or its prediction as the only one there is.

Returns (n_data, n_realizations) either way, so that whatever reads it does not have to care which it got.

geoml.plots.prepare.realization_store(var)

A variable’s realizations as a store to read in bands, not as an array.

The same (n_data, n_realizations) that realizations hands back, without asking memory for all of it at once. A block model’s simulations are the one thing a container holds that will not fit – hundreds of gigabytes is an ordinary size for them – and every use of them here is a reduction over locations, which never needs more than a band of rows at a time.

A variable with no simulations falls back on its prediction, which is a single column and already in RAM, and comes back as a one-band store so that the caller has only the one path to write.

geoml.plots.prepare.block_density(container, density, n_data, n_realizations)

A density per block, to turn a volume into a tonnage.

density may be a number, the name of a metadata column, or the name of a ContinuousVariable. Only the last of these can be uncertain, and when it is, its realizations are matched one to one with the grade’s: simulation i of the density belongs with simulation i of the grade, and pairing them any other way would invent a correlation that was never modelled.

Comes back as a single number where the density is one, and otherwise as a store to read in bands alongside the grade. A simulated density is exactly as big as the grade, so materializing it would give back everything not materializing the grade saves.

geoml.plots.prepare.containing_variable(container, var)

The variable var is a component of, if it is a component of one.

geoml.plots.prepare.uncertainty_values(container, name, variables=())

A number per location saying how much the model doubts itself there.

Which number that is depends on what was modelled, so it is named rather than guessed at. name may be:

• an array, one value per location. Whatever the number is and wherever it came from, it can always be handed over directly, which is the way out of every case the names below do not reach.

• a path, naming exactly where to read it – “Elements/uncertainty” while the grade is one of its components, or a column belonging to some other variable entirely. (The old dotted “Variable.column” is refused with the replacement spelled out: a . inside a label would make a wrong guess look like a working one.)

• a bare name, looked for on each of variables in turn and then among the metadata columns.

variables is what a bare name is tried against: the grade, and then the variable containing it. That second one is not a nicety. A component of a vector variable has latent_variance set to None and no uncertainty of its own – the column that exists belongs to the parent, so grading Zn and asking for “uncertainty” has to reach Elements to find anything.

geoml.plots.prepare.block_volume(container)

What one block is worth, as a number or as one value per block.

A regular grid has a single spacing and so a single volume; a BlockSet3D carries a size per block and answers with a column. Kept apart from the density because only one of them is a property of the container.

geoml.plots.prepare.grade_tonnage(container, name, density=None, cutoffs=30, uncertainty=None, max_uncertainty=None)

How much material sits above a cut-off, and how good it is.

Each block contributes its volume, or its mass where a density is given, to every cut-off its grade clears. Simulations are carried through separately rather than averaged first: the curve of the mean model is not the mean of the curves, since a cut-off is a threshold and averaging either side of it gives different answers.

The simulations are read a band of blocks at a time and never held whole: a block model runs to hundreds of gigabytes of them, and what comes out is one small number per cut-off per realization. Each block is placed at the highest cut-off it clears and the curve is the running total from the top down, so the cost is one pass over the grade rather than one per cut-off. Giving cutoffs as values rather than as a count saves the pass that would otherwise be needed to find their range.

Parameters:

• container – A gridded container – the volume of a block comes from its spacing.

• name (str) – The variable to take as the grade.

• density (float or str) – A number, a metadata column, or a ContinuousVariable. Without one the curve is in volume.

• cutoffs (int or array-like) – The grades to cut at, or how many of them to spread evenly across the range of the data.

• uncertainty (str or array) – Where to read how sure the model is at each block: a column name, a a path (“Variable/column”) naming which variable it belongs to, or the values themselves. See uncertainty_values. A bare name is looked for on the grade, then on the variable containing it, then in the metadata.

• max_uncertainty (float) – Blocks doubted more than this are left out altogether – not counted at any cut-off, and not counted towards the grade above one. A block the model cannot speak for is not tonnage.

Returns:

dict – cutoff (n_cutoffs,); tonnage, grade and metal (n_cutoffs, n_realizations); unit, which is the extent of a block or “mass” depending on whether a density was given; and kept and total, the blocks that survived the uncertainty filter and the blocks there were.

geoml.plots.prepare.grid_shape(n_panels)

Rows and columns for n_panels, as square as they go.

Style

One place defining what a geoML figure looks like.

Every function in this sub-package draws inside context(), a matplotlib.rc_context, so the settings below apply to the figure being drawn and to nothing else: importing geoML does not change how anyone else’s plots look. Call use() to take them globally, if that is what you want.

PALETTE is what tells one category from another, and it is the same list wherever a category is drawn, so a rock type keeps its colour from a histogram to a map to the dashboard. Change it here and every figure follows.

geoml.plots.style.context()

The geoML settings, for the figure being drawn and nothing else.

geoml.plots.style.use()

Takes the geoML settings globally, for every figure from here on.

geoml.plots.style.color(index)

The colour of the index-th category, cycling when there are many.