Okay, things are a bit clearer with your updated question. Thoughts that occur:
- How dense is the data within each cluster compared to the rest of the space?
In your diagrams the square regions all look separated by fairly wide gaps. If you have enough samples within each cluster, and you have few/zero outliers, then you might be able to just ignore the shape prior and cluster based on pairwise sample distances. A graph cut type approach, or spectral clustering, for example.
- Have you tried simple clustering algorithms but with different distance metrics?
You rectangular regions all look axis aligned, yet you talk about euclidean distance clustering. Have you tried tweaking an algorithm to use $L_1$ distances instead? It won't be ideal (some of the above shapes are cuboids, not cubes), but if your regions are suitable spread out, it might be good enough.
- Have you tried parametric modelling?
It's going to be fairly easy to define a parametric model for a cuboid, and an associated likelihood for samples based on this (e.g. 1/(vol. of cuboid) for samples within it, zero otherwise). You could do this for a mixture model of cuboids. From here you could take an optimisation route (throw the log likelihood, plus any log prior or regularisation terms you desire, into your favourite optimisation algorithm) or you could do the whole bayesian model selection thing, marginalising over the parameters. My biggest worry here is that if your dimensionality is high, you may end up with a lot of parameters to deal with.
Any of that sound doable, or are there issues with it I've missed?
EDIT: As requested, more details on parametric modelling.
Wikipedia isn't a bad place to start (http://en.wikipedia.org/wiki/Parametric_model).
Let's talk about this particular problem though. It's easiest to think about it in 1D at first. Let's say we have a uniform distribution in a 1d space, specified by a start point $a$ and an end point $b$. These are our parameters. If we observe a sample $x$ in our 1D space, the probability of that sample given $a$ and $b$ is
$$
p(x|a,b) = \begin{cases} \frac{1}{|b-a|} & \text{if } a \leq x < b \\ 0 & \text{ otherwise}\end{cases}
$$
which is also known as the likelihood of the parameters.
We can extend this beyond one dimension very easily. In that case, $a$ and $b$ become opposite corner vertices in our space, the test $a \leq x < b$ is performed in each dimension separately, and $|b-a|$ become the volume of the cuboid they describe.
E.G. for a 2D space maybe $a = (2,2)$ and $b = (4,5)$. Any $x$ which lies within the cuboid they specify would have a $p(x|a,b)$ equal to $1/((4-2)*(5-2)) = 1/6$. A diagram is probably makes this clearer...
_____________.(4,5)
| |
| p(x|a,b) = |
| 1/6 | p(x|a,b) = 0
(2,2).|____________|
You, however, don't have one cuboid cluster. You have a number of them. So it makes sense to use a mixture model (https://en.wikipedia.org/wiki/Mixture_model). A mixture model is where we take a group of simple distirbutions and add them up to make a more complicated one - just like the way the distribution you're interested in is equal to the sum of lots of cuboid distributions. Let's start with a two component model. In this situation you have two cuboid clusters in your space, so we need to specify four corner parameters $a_1, b_1, a_2, b_2$ (two for each cuboid) and two weighting parameter $w_1$ and $w_2$ which are always positive and sum to 1 (to ensure our resulting distribution will integrate to 1). Then we have
$$
p(x|a_1, b_1, a_2, b_2, w) = w_1*p(x|a_1,b_1) + w_2*p(x|a_2,b_2)
$$
where $p(x|a_1,b_1)$ is just the likelihood for a single cuboid we defined earlier. The $w$ parameters control how common it is for a sample to be part of a particular cluster, and the $a$ and $b$ parameters control where that particular cluster is. If we add another cluster to our earlier example and set $w_1 = w_2 = 0.5$ we'd have something which looks like this:
______________.(4,5)
| |
|w1*p(x|a1,b1)| p(x|a1,b1, w1, a2, b2, w2) = 0
| =1/12 |
(2,2).|_____________| ________________________________. (2,10)
| w2*p(x|a2,b2) = 1/8 |
(1,6).|________________________________|
which is already starting to look a bit like your examples.
- It's easy to extend this from two components to $K$ components:
$$
p(x|a_1, b_1, w_1 \ldots a_K, b_K, w_K) = \sum_{k=1}^K w_k*p(x|a_k, b_k)
$$
where $w_k\geq 0$ and $\sum_{k=1}^K w_k = 1$.
- It's also easy to extend this to many observed sample instead of just 1. If you have $N$ samples, then assuming they are iid (https://en.wikipedia.org/wiki/Iid) given the parameters, we have
$$
p(x_1 \ldots x_N|a_1,b_1, w_1, \ldots a_K, b_K, w_K) = \prod_{i=1}^N p(x_i|a_1,b_1, w_1, \ldots a_K, b_K, w_K).
$$
That, then, gives you the likelihood of a $K$ component cuboid mixture model, evaluated on a particular dataset of $N$ samples. You could try optimising that directly to get the maximum likelihood result (protip - optimise the log likelihood, and include some regularisation to stop the cuboids getting arbitrarily small). You could also include a prior over the parameters to get the MAP result.
You may run into optimisation problems due to the discontinuities at the edge of each cuboid. You're going to need to model outliers somehow too,if they occur. You're also going to have to find a way to select the number of clusters, $K$. And actual performance is going to be heavily dependant on how you initialise your parameters. But these are problems to think about once you've decided if this is actually a route you want to go down (as opposed to something more out of the box), and have checked how much of a problem they cause in practise.
The only reference I can find right now to prior work (though I'm sure people have done this many times before):
Parameter estimation for finite mixtures of uniform distributions
Craigmile, PF (Craigmile, PF); Titterington, DM (Titterington, DM)
Source: COMMUNICATIONS IN STATISTICS-THEORY AND METHODS Volume: 26 Issue: 8 Pages: 1981-1995 DOI: 10.1080/03610929708832026 Published: 1997
Also that took way more time to explain than I anticipated, which makes me wonder if writing the code might be too big an investment of your time too if it's all new to you :P. Try R-trees first, as other people have suggested.