2017-05-29 4 views
2

Ich mache Simulationsstudien am linearen Modell y = X * beta + eps mit der Größe (X) = [n d]. Ich betrachte den Effekt der Dimensionalität d basierend auf zwei Methoden. Ich führe 10 simulierte Daten und bekomme entsprechende Beta-Schätzungen, und dann Ich möchte den Mittelwert der Beta über die 10 simulierten Daten berechnen.Zell-Array-Mittelwert in Matlab

Mein MATLAB-Code Spielzeug ist wie folgt:

 nsim=10; %iteration number 
     dd=[4 6]; %two dimension cases,\beta=(\beta_1,\cdots,\beta_d)^T 
     ddlen=length(dd); 
     nmethod=2; %two methods 
     seednum=0; 

     BH = cell(nsim,ddlen,nmethod); %estimation of beta w.r.t two dimension cases and two methods 

     for di = 1:ddlen 
      d = dd(di); 
      for ni = 1:nsim 
       seednum = seednum + di*ni; 
       randn('seed', seednum); 
       betahat=randn(d,1); 
       for method = 1:nmethod 
        if method==1 
         BH{ni,di,method} = betahat; 
        else 
         BH{ni,di,method} = 10*betahat; 
        end 
       end 
      end 
     end 

Dann können wir

BH(:,:,1) = 

    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 


BH(:,:,2) = 

    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 
    [4x1 double] [6x1 double] 

erhalten möchte ich den Mittelwert über 10 Zeilen (NSIM = 10) und erhalten so etwas wie

mean(BH(:,:,1))= 

    [4x1 double] [6x1 double] 

mean(BH(:,:,2)) = 

    [4x1 double] [6x1 double] 

Irgendwelche Ideen? Vielen Dank!

+0

Thanks @ EBH. Aber deine Antwort ist nicht, dass ich will. Die Rückkehr sollte zwei Vektoren sein, einer ist [4x1 double] und der andere ist [6x1 double] und mit anderen Worten, der Mittelwert von 10 [4x1 double] bzw. der Mittelwert von 10 [6x1 double]. –

+3

Warum verwenden Sie Zellenarrays, wenn alle Ihre Arrays die gleiche Größe haben? – beaker

+0

Warum legen Sie einen neuen Zufalls-Seed in jeder Schleife fest? – EBH

Antwort

0

Ich weiß nicht, ob es die effizienteste Methode ist, dass Sie tun, aber Sie arrayfun dafür verwenden können:

% generate random array 
BH = repmat({rand(4,1),rand(6,1)},[10 1 2]); 
% generate indexes for the 2nd and 3rd dimensions 
[n2,n1] = meshgrid(1:size(BH,2),1:size(BH,3)); 
% get mean across 1st (cell) dimension 
[res] = arrayfun(@(n1,n2)mean([BH{:,n1,n2}],2),n1(:),n2(:),'UniformOutput',false); 
% reshape to desired output 
res = reshape(res,[1 size(BH,2) size(BH,3)]); 

Falls Sie zu einem N-dimensionalen Zellenfeld verallgemeinern wollen:

% generate random array 
BH = repmat({rand(4,1),rand(6,1)},[10,1,2,2,5]); 
sz = size(BH); 
% generate indexes for the 2nd and 3rd dimensions 
n = cell(1,numel(sz) - 1); 
[n{:}] = ndgrid(1:sz(2),1:sz(3),1:sz(4),1:sz(5)); 
n = cell2mat(cellfun(@(x) {x(:)},n)); 
idx = 1:size(n,1); 
% get mean across 1st (cell) dimension 
[res] = arrayfun(@(idx)mean([BH{:,n(idx,1),n(idx,2),n(idx,3),n(idx,4)}],2),... 
    idx,'UniformOutput',false); 
% reshape to desired output 
res = reshape(res,[1 sz(2:end)]); 
+0

Es funktioniert! Danke @ user2999345. Nun, ich finde es ein bisschen kompliziert, eine Verallgemeinerung zu machen, die Ihrem Ansatz entspricht. Denken Sie nur an BH = repmat ({rand (4,1), rand (6,1)}, [10,1,2,2,5]); dh 10 Iterationen, 2 Dimensionen, 2 Geräuschpegel und 5 Kandidaten Methoden für lineares Modell y = X * beta + eps, wobei eps ~ N (0, \ sigma^2 * I) und \ sigma ist der Geräuschpegel –

+0

check out my edit, ich hoffe, es beantwortet Ihre Anfrage – user2999345

+0

Wow , deine neue Bearbeitungsantwort ist fantastisch! Danke vielmals! @ user2999345 –

0

Wenn ich dich richtig verstehe, willst du den Mittelwert über alle Elemente nehmen, die sich an derselben Position im Vektor befinden. So erhalten wir aus allen Vektoren in BH(:,1,1) einen Vektor von 4 Mitteln, jeweils für eine Position im Vektor. Gleiches gilt für BH(:,1,2). Für BH(:,2,1) und BH(:,2,1) machen wir dasselbe, aber mit 6 Elementen in einem Vektor.

Sie können den folgenden Code für das verwenden:

% split BH to 2 arrays: 
bh4 = reshape(cell2mat(BH(:,1,:)),[],nsim,2); % all the 4 elements vectors 
bh6 = reshape(cell2mat(BH(:,2,:)),[],nsim,2); % all the 6 elements vectors 
meanBH4 = squeeze(mean(bh4,2)); % mean over all 4 element vectors 
meanBH6 = squeeze(mean(bh6,2)); % mean over all 6 element vectors 

jedoch ein Schritt in dem richtigen Weg, dies zu tun wäre, zwei Arrays zu definieren, eine für jede Methode:

BH1 = zeros(nsim,ddlen,dd(1)); 
BH2 = zeros(nsim,ddlen,dd(2)); 

dann in Ihrer Schleife die Werte zu ihnen zuweisen:

if method==1 
    BH1(ni,di,:) = betahat; 
else 
    BH2(ni,di,:) = 10*betahat; 
end 

und am Ende nehmen Sie einfach den Mittelwert von jedem:

meanBH1 = mean(BH1,3) 
meanBH2 = mean(BH1,3) 

EDIT:

Um all dies in einer 'Matlabish' Art und Weise schreibe ich folgendes vorschlagen würde:

nsim = 10; % iteration number 
dd = [4 6]; % two dimension cases,\beta=(\beta_1,\cdots,\beta_d)^T 
methods = 2; % two methods 

% preapering random seeds 
s = bsxfun(@times,1:numel(dd),(1:nsim).'); 
seednum = cumsum(s(:)); 

% initialize results array 
BH = nan(max(dd),nsim,numel(dd),methods); 
counter = 1; 
for k = 1:numel(dd) 
    for n = 1:nsim 
     % set a new random seed from the list: 
     rng(seednum(counter)); 
     % compute all betahats with this seed: 
     betahat = randn(max(dd),2).*repmat([1 10],[max(dd) 1]); 
     % assign the values to BH by dimension: 
     for m = 1:methods 
      BH(1:dd(k),n,k,m) = betahat(1:dd(k),m); 
     end 
     counter = counter+1; 
    end 
end 
% compute the means over iterations: 
means = squeeze(mean(BH,2,'omitnan')) 

und so erhalten Sie means als Ihr Ergebnis .


P.S. Ich weiß nicht, warum Sie randn('seed', seednum) bei jeder Iteration nennen, besides that's not a recommended syntax, aber Sie, wenn Sie löschen können, dann können Sie die meisten der Schleifen vektorisiert, und Ihr Code quetscht auf:

% compute all betahats: 
betahat = randn(nsim,max(dd),numel(dd),2); 
% apply dimensions: 
for k = dd 
    betahat(:,k+1:end,1,:) = nan; 
end 
% apply methos 2: 
betahat(:,:,:,2) = betahat(:,:,:,2)*10; 

% compute the means over iterations: 
means = squeeze(mean(betahat,1,'omitnan')) 

Hoffnung Dinge jetzt klarer sieht. ..

+0

Danke @ EBH. Aber deine Antwort ist nicht, dass ich will. Die Rückkehr sollte zwei Vektoren sein, einer ist [4x1 double] und der andere ist [6x1 double] und mit anderen Worten, der Mittelwert von 10 [4x1 double] bzw. der Mittelwert von 10 [6x1 double]. –

+0

@JohnStone, was Sie in der Frage beschreiben, ist eine 1 * 2 * 2 Ausgabe, und nicht zu Vektoren von 4 * 1 und 6 * 1 Ausgabe. Siehe jedoch meine Bearbeitung dafür. – EBH

+0

Danke @ EBH. Nun, wenn ich eine komplexere Einstellung betrachte, wie etwa BH = repmat ({rand (4,1), rand (6,1)}, [10,1,2,2,5]); dh 10 Iterationen , 2 Dimensionen, 2 Rauschpegel und 5 Kandidatenmethoden für das lineare Modell y = X * beta + eps, wobei eps ~ N (0, \ sigma^2 * I) und \ sigma der Rauschpegel ist. Was sollte ich tun? –

0

, Alternativ

% split into seperate cell arrays 
BH_1 = BH(:,:,1); 
BH_2 = BH(:,:,2); 

% create matrix of compatible vectors, and take mean and put result back into cell array 
BH_1_mean = cat(2,{mean(cell2mat(BH_1(:,1)'),2)}, {mean(cell2mat(BH_1(:,2)'),2)}); 
BH_2_mean = cat(2,{mean(cell2mat(BH_2(:,1)'),2)}, {mean(cell2mat(BH_2(:,2)'),2)}); 
+0

Danke @ kedarps, es funktioniert! –