The Voronoi texture uses a square grid, which is apparent at lower values of randomness. However, is it possible to use a hexagonal grid instead? Adding the right linear transformation could make the points' locations align with a hexagonal grid, but the circles around each point would distort into ellipses in doing so. More precisely, I would like to reproduce the below image (not made in Blender).
A Voronoi texture works by dividing the texture-space into cells, defining a pseudo-random 'feature point' in each cell, which is a function of the cell's location. Then each shading-point looks up the the feature-point in its own and neighbouring cells, returning the distance to the closest feature-point.
Blender's Voronoi works on the basis of rectangular cells. The ask, here, is to base a Voronoi on hexagonal cells.
1. Create a hexagonal grid
That's dealt with at the bottom of this answer. One of the outputs is '2D Index' which is the coordinate in the current texture-space of the centre of each cell.
2. Create a feature-point in each cell
To obtain random points, evenly distributed in a hexagon, the hexagon is divided into 3 rhombuses as shown below. One of those is randomly selected using the X of 3D White Noise. Then random 0-1 multiples,(using Y and Z of the noise) of the U and V vectors for that rhombus yield evenly distributed points inside it. The result is shown on the right:
The nodes for this:
A driven 'Randomness' mix is included to blend between cell-centers and random offsets of them. This should be wired in, but none of this tree is optimized.
3. For each cell, return the coordinates of its own, and its neighbour's centers
(above). '0' is the cell itself, '1' is the cell to the right, and the others are obtained by rotating '1' through multiples of 60 degrees. Again, to optimise, the rotations should be avoided: the resulting vectors should be hard-coded. In this version, I'm avoiding spaghetti at the cost of efficiency.
4. For each shading-point, find the nearest feature point
A small helper-group returns the closer of the points A and B to P:
which is daisy-chained to ...
5. Yield the coordinate of the closest feature-point to the shading point
... when used with the combinations already described.
None of this is optimised for speed. Rotations are used instead of hard-coding, distances are repeatedly calculated instead of being passed on.. and so on. The White Noise is very sensitive to floating-point errors in input values around -1, it seems? So a large number has been added to the input vectors, when it's used, here.
This .gif shows 'Distance', 'Position', and 'Color' outputs of the whole kaboodle, with 'Randomness' being swiped from 0-1:
Settings for your ref: about .7 'Randomness', and a color-ramp:
Voronoi noise is made by tiling circles in a square fashion and randomising the positions of all circles. Since the circles are limited to their tile's size, when you randomly change the position, they exceed their tile size and start to form seams. To eliminate that, we check the neighboring tiles i.e left, right, up, down, etc. This is how the usual voronoi is made. In the case of hexagonal voronoi, we only need to change one thing. We need to tile circles using hexagonal tiles / hexagonal coordinates. You can see how to make them here How To Do UV Indexing in hexagonal pattern? OR you can use the "hextiles" nodegroup in the blend file below. I have made a video how I made this : https://www.youtube.com/watch?v=RbH3_d9j7ro
Here is the blend: