2016-02-18 1 views
2

Kind all,Représentation graphique de la densité des itérations

Je travaille dans MATLAB et j'utilise les techniques de Monte Carlo pour ajuster un modèle. En fait, si nous supposons que mon modèle est une fonction simple comme

y = m * x^2 + c

Et que mes deux paramètres m et c varient entre 0,5 et 10, je peut tirer aléatoirement d'un tel espace de paramètre et obtenir à chaque fois un nouveau y. Si je conspire toutes mes réalisations de y-je obtenir quelque chose de similaire à la figure suivante:

Multiple Iterations of y, colours are drawn by the hsv colormap

Y at-il une façon de représenter la densité des réalisations? Je veux dire, y a-t-il un moyen (au lieu de tracer toutes les réalisations) d'obtenir une sorte de tracé de contour qui se situe entre le minimum de mes itérations et le maximum pour lequel sa couleur représente la quantité de réalisations comprises dans un certain intervalle?

Merci à tous!

+1

Comment définissez-vous la densité? Vous avez besoin d'une sorte de définition de cela pour être en mesure de réaliser ce que vous voulez –

Répondre

0

Vous pouvez calculer y pour des points discrets de x, tout en définissant des valeurs aléatoires pour c et m. Ensuite, en utilisant la fonction hist, vous pouvez trouver une "densité non normalisée" de valeurs de fonction pour un x donné. Vous pouvez ensuite le normaliser pour obtenir une densité réelle des valeurs, de sorte que la zone globale sous la courbe de distribution résume jusqu'à 1.

Pour le visualiser, vous construire une grille de maille [X, Y] le long des valeurs de x et y et mettre les valeurs de densité comme Z. Maintenant, vous pouvez soit tracer le surf ou son tracé de contour.

enter image description here

enter image description here

enter image description here

Voici le code:

clear; 

n = 1000000; %number of simulation steps 

%parameter ranges 
m_min = 0.5; m_max = 10; 
c_min = 0.5; c_max = 10; 

%x points 
x_min = 1; x_max = 4; x_count = 100; 
x = linspace(x_min, x_max, x_count); 
x2 = x.^2; 

y_min = 0; y_max = m_max*x_max*x_max + c_max; y_step = 1; 

m = rand(n, 1)*(m_max - m_min) + m_min; 
c = rand(n, 1)*(c_max - c_min) + c_min; 

c = repmat(c, 1, x_count); 

y = m*x2 + c; 

x_step = (x_max- x_min)/(x_count-1); 
[X, Y] = meshgrid(x_min:x_step:x_max, y_min-y_step:y_step:y_max+y_step); 
Z = zeros(size(X)); 

bins = y_min:y_step:y_max; 

for i=1:x_count 

    [n_hist,y_hist] = hist(y(:, i), bins); 

    %add zeros on both sides to close the profile 
    n_hist = [0 n_hist 0]; 
    y_hist = [y_min-y_step y_hist y_max+y_step]; 

    %normalization 
    S = trapz(y_hist,n_hist); %area under the bow 
    n_hist = n_hist/S; %scaling of the bow 

    Z(:, i) = n_hist'; 
end 

surf(X, Y, Z, 'EdgeColor','none'); 
colormap jet; 
xlim([x_min x_max]); 
ylim([y_min y_max]); 
xlabel('X'); 
ylabel('Y'); 

figure 
contour(X,Y,Z, 20); 
colormap jet; 
colorbar; 
grid on; 
title('Density as function of X'); 
xlabel('X'); 
ylabel('Y'); 

Une autre vue intéressante est un tracé de chaque section en fonction de la valeur x:

enter image description here

Voici le code pour ce terrain:

clear; 

n = 1000000; %number of simulation steps 

%parameter ranges 
m_min = 0.5; m_max = 10; 
c_min = 0.5; c_max = 10; 

%x points 
x_min = 1; x_max = 4; x_count = 12; 
x = linspace(x_min, x_max, x_count); 
x2 = x.^2; 

m = rand(n, 1)*(m_max - m_min) + m_min; 
c = rand(n, 1)*(c_max - c_min) + c_min; 

c = repmat(c, 1, x_count); 

y = m*x2 + c; 

%colors for the plot 
colors = ... 
[ 0 0 1; 0 1 0; 1 0 0; 0 1 1; 1 0 1; 0 0.75 0.75; 0 0.5 0; 0.75 0.75 0; ... 
1 0.50 0.25; 0.75 0 0.75; 0.7 0.7 0.7; 0.8 0.7 0.6; 0.6 0.5 0.4; 1 1 0; 0 0 0 ]; 

%container for legend entries 
legend_list = cell(1, x_count); 

for i=1:x_count 
    bin_number = 30; %number of histogramm bins 
    [n_hist,y_hist] = hist(y(:, i), bin_number); 
    n_hist(1) = 0; n_hist(end) = 0; %set first and last values to zero 

    %normalization 
    S = trapz(y_hist,n_hist); %area under the bow 
    n_hist = n_hist/S; %scaling of the bow 
    plot(y_hist,n_hist, 'Color', colors(i, :), 'LineWidth', 2); 
    hold on; 

    legend_list{i} = sprintf('Plot of x = %2.2f', x(i)); 
end 

xlabel('y'); 
ylabel('pdf(y)'); 
legend(legend_list); 
title('Density depending on x'); 
grid on; 
hold off; 
1

Ce n'est pas très joli, mais vous pouvez faire varier les paramètres et jouer avec le scatter/traçage, pour le rendre un peu plus attrayant. J'ai aussi supposé une distribution gaussienne au lieu de votre distribution aléatoire (une coloration totalement aléatoire vous donnera une densité uniforme). Aussi ce code pourrait être optimisé pour la vitesse.

n = 1000; 
l = 100; 
x = linspace(1, 10, l); 
y = repmat(x.^2, n, 1); 
c = repmat(normrnd(1, 1, n, 1), 1, l); 
m = repmat(normrnd(1, 1, n, 1), 1, l); 

y = y.*m + c; 
p = plot(y', '.'); 
figure; hold on; 
for i = 1:l 
    [N,edges,bin] = histcounts(y(:, i)); 
    density = N./sum(N); 
    c = zeros(n, 3); 
    for j = 1:n 
     c(j, :) = [1-density(bin(j))/max(density), 1-density(bin(j))/max(density), 1-density(bin(j))/max(density)]; 
    end 
    scatter(ones(n, 1)*x(i),y(:, i),[],c,'filled'); 
end 

Confère

enter image description here

Cela crée un histogramme des valeurs de y pour chaque position x, calcule alors la densité de probabilité pour chaque valeur de y et de couleurs dans les points. Ici, les valeurs y de chaque position x sont normalisées individuellement, pour colorer les points en fonction de la densité globale de l'intrigue dont vous avez besoin pour renormaliser.