Generate a vector under norm condition

20 views (last 30 days)
Youssef Fakih
Youssef Fakih on 22 May 2020
Commented: Rik on 2 Mar 2021
I want to Generate a random vector :
A = randn(10,1);
Such as the Norm of A is less than a constant value "For example 0.01"
norm(A) < 0.01
How to Generate this vector under a norm condition ?
Do I have to use while loop to test if the condition is satisfied every time ?
Or there is some predefined function that can perform such type of condition ?
Thank you for your help
  2 Comments
Rik
Rik on 22 May 2020
Are you sure you actually mean the vector norm? That seems strange in combination with randn. Since the norm function will return the square root of the sum of the squares, a value as low as 0.01 seems fairly strict.
norm_val=zeros(1,10^6);%run a million trials
for n=1:numel(norm_val)
A = randn(10,1);
norm_val(n)=norm(A);
end
figure(1),clf(1)
histogram(norm_val)
title(sprintf('the norm was below 0.01 %d times with %.0e trials',sum(norm_val<0.01),numel(norm_val)))
John D'Errico
John D'Errico on 22 May 2020
Edited: John D'Errico on 22 May 2020
As Rik points out, this is a TINY ball for the sample to live in.
We can use the chi square distibution to see the probability of such an event happening just by random chance. A chi-square distribution is the distribution of the sum of squares of N normally distributed random variables. As it turns out, this gives us what we need to compute the probability of such an event occuring randomly, as
chi2cdf(0.01^2,10)
ans =
2.6041e-24
I used 0.01^ there, because you want the norm to be that small. If the norm of the vector is that small, then the sum of the squares of the elements must be less than 0.01^2.
We see a REALLY tiny chance of that happening randomly. All of the elements of the final vector will be almost identcally zero. You REALLY don't want to use rejection sampling on something like this.
As such, I wonder if this is the real question, but my answer will show how to solve it, for some arbitrary maximum value of Euclidean norm with a normal vector.

Sign in to comment.

Answers (2)

John D'Errico
John D'Errico on 22 May 2020
Edited: John D'Errico on 22 May 2020
As long as you are willing to specify that the vector is normally distributed under a constraint, you can do this. In fact, it is even pretty easy. One point you did not mention is if this is the Eiclidean norm. That is usually implied when all that is stated is the word norm though.
You are essentially asking for what might be called a truncated multi-variate normal distribution. Here, the truncation is such that the point in a 10-dimensional space must live inside a sphere of radious 0.01. Any point inside that sphere is acceptable. Now, admittedly, a sphere of radius 0.01 in a 10 dimensional space is TINY. Almost every element of that vector wil be essentially zero.
As I said in my comment to your question, we can think of this as selecting a random point based on a truncated normal, such that the solution lies inside a unit ball of radius R. In your question, you used R = 0.01.
I'll make the problem a little more general, and solve it for vectors of length n. Again, in your question, you had n=10. The final question is how many of these vectors do you wish to generate? You said you only want 1 vector, but I'll assume you want to generate nvec such vectors.
function X = normalballsample(n,R,nvec)
% generates a truncated normal array of size (n,nvec) such that the 2-norm of the columns does not exceed R
% if nvec was not provided, just generate one vector.
if nargin < 3
nvec = 1;
end
% generate a uniform sample, of size nvec, that lies between 0 and chi2cdf(n,r^2)
u = rand(1,nvec)*chi2cdf(R^2,n);
% invert those numbers through the chi-square CDF
% then take the sqrt. This uses the classic inverse CDF method to generate numbers
% from the desired distribution. Take the sqrt, since we really wanted a chi sample.
r = sqrt(chi2inv(u,n));
% generate an unconstrained sample in n dimensions
X = randn(n,nvec);
% scale the samples so they ALL have exactly unit norm, then rescale them
% so the norm will be r.
xnorms = sqrt(sum(X.^2,1));
% in case someone has an old MATLAB release than R2016b, use bsxfun
% X = bsxfun(@times,X,r./xnorms);
X = X.*(r./xnorms);
end
Does this work? Well, it should. ;-) It will take some testing to convince you, perhaps.
First, generate a sample in 2 dimensions, such that all points lie inside a norm of 3. Generate 10000 such vectors of length 2, and plot.
X = normalballsample(2,2,10000);
plot(X(1,:),X(2,:),'.')
axis equal
grid on
xlabel X
ylabel Y
So the points must lie inside the unit circle, centered at the origin. They look normally distributed, with a truncation limit at the perimeter of a circle of radius 2.
How does it work on a vector of length 10, constrained to lie inside a norm of 0.01?
X = normalballsample(10,.01,1)
X =
-0.00129143216137951
-4.59296749685007e-05
0.000972209117713196
0.00608141153128656
-0.00287329676048124
0.00296727125014177
-0.00348917693949026
0.00318993726271893
-0.00154208722267966
-0.00384775092173714
norm(X)
ans =
0.00980777609363417
As expected, the norm will usually be pretty close to 0.01 when the radius is that small.
The above code uses tools from the stats toolbox, but it could easily enough have been written to avoid that toolbox, in case that toolbox is not available. The chi-squared distribution should be convertable via transformation to a gamma distribution. But I hate to waste brain sweat if it can be avoided.
As well, the code I wrote assumes R2016b or later. But the one line that would be a problem is indicated, and I show how to write is using bsxfun.
  3 Comments
Lukas Voss
Lukas Voss on 2 Mar 2021
I am trying to solve a similar problem. My goal is to generate nvec random vectors that are located in a spherical shell around the origin, so between some r_min and r_max. How should I modify your function to get this done?
Thanks in advance!
Rik
Rik on 2 Mar 2021
You can either do something clever, or just take this function, let it generate values, and remove all points that have a norm below r_min.

Sign in to comment.


Abdolkarim Mohammadi
Abdolkarim Mohammadi on 22 May 2020
I don't know of any nice functions to do that in one shot, but I think you can start with a random point and use a solver to optimize it for you so that it brings the norm to the desired value. But be aware that the data may get non-Gaussian.
  1 Comment
John D'Errico
John D'Errico on 22 May 2020
The sample will then lie ON the unit ball. So the norm would then be 0.01, not less than or equal to 0.01.

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!