8

我有一个 5000x5000 的网格,我正在尝试在 MATLAB 中实现一个简单的癌症划分模型。最初,它选择一个随机点 (x,y) 并使该细胞成为癌细胞。在第一次迭代中,它分裂 - 父细胞留在原地,子细胞被随机分配给任何相邻细胞。
到目前为止很容易。

我的问题是这样的:在连续的迭代中,一个子细胞通常会被分配给一个已经有癌细胞的细胞。在这种情况下,我希望子单元格取代它的位置并将已经存在的单元格“撞击”到相邻的单元格。如果该相邻单元格为空,则将其填充并且该过程停止。如果没有,则已经在该位置的单元格被碰撞,依此类推,直到最后一个单元格找到一个空白空间并且该过程停止。

这应该很简单,但我不知道如何编写代码以及使用什么样的循环。
我是物理科学家而不是程序员,所以请把我当傻子一样对待!

4

2 回答 2

4

这是我一起破解的一个功能,大致符合您提供的规格。

随着癌细胞的数量变大,我确实放慢了速度。

基本上我有几个变量,表示单元格位置网格的 NxN 矩阵(我称之为 a plate,因为网格是现有 matlab 函数的名称)

我可以快速迭代的点向量。我选择一个种子位置,然后运行一个 while 循环,直到网格已满。

在每次循环迭代中,我对每个单元格执行以下操作:

  • 生成一个随机数以确定该单元格是否应该分裂
  • 生成随机方向进行划分
  • 找到该方向的第一个开盘位置
  • 填充该位置

我没有对它进行广泛的测试,但它似乎有效。

function simulateCancer(plateSize, pDivide)

plate = zeros(plateSize, plateSize);
nCells = 1;
cellLocations = zeros(plateSize*plateSize,2);

initX = randi(plateSize);
initY = randi(plateSize);

cellLocations(nCells,:) = [initX, initY];

plate(initX, initY) = 1;

f = figure;
a = axes('Parent', f);
im = imagesc(plate, 'Parent', a);


while(nCells < (plateSize * plateSize))
    currentGeneration = currentGeneration+1;
    for i = 1:nCells
        divide = rand();
        if divide <= pDivide
            divideLocation = cellLocations(i,:);
            divideDir = randi(4);
            [x, y, v] = findNewLocation(divideLocation(1), divideLocation(2), plate, divideDir);
            if (v==1)
                nCells = nCells+1;
                plate(x,y) = 1;
                cellLocations(nCells,:) = [x,y];
            end
        end
    end
    set(im,'CData', plate);
    pause(.1);
end

end

function [x,y, valid] = findNewLocation(xin, yin, plate, direction)   
    x = xin;
    y = yin;
    valid = 1;
    % keep looking for new spot if current spot is occupied
    while( plate(x, y) == 1)
       switch direction
            case 1 % divide up
                y = y-1;
            case 2 % divide down
                y = y+1;
            case 3 % divide left
                x = x-1;
            case 4 % divide down
                x = x+1;
            otherwise
            warning('Invalid direction')
            x = xin;
            y = yin;
        return;
       end

       %if there has been a collision with a wall then just quit
       if y==0 || y==size(plate,2)+1 || x==0 || x==size(plate,1)+1 % hit the top
           x = xin; %return original values to say no division happend
           y = yin;
           valid = 0;
           return;
       end

    end


end

注意:我没有考虑推送单元格,而是以一种将单元格留在当前位置并在行/列末尾创建新单元格的方式对其进行了编码。只要您不关心世代,它在语义上不同但在逻辑上具有相同的最终结果。

于 2012-06-22T14:38:33.710 回答
2

另一个问题的启发,我想使用图像处理技术来实现这个模拟。具体来说,我们可以使用形态扩张来扩散癌细胞。

这个想法是使用如下结构元素来扩展每个像素:

1 0 0
0 1 0
0 0 0

其中中心是固定的,另一个1随机放置在其余八个位置之一。这将有效地沿该方向扩展像素。

执行膨胀的方式是创建一个空白图像,只有一个像素集,然后使用简单的 OR 运算累积所有结果。

为了加快速度,我们不需要考虑每个像素,只需要考虑由癌细胞簇形成的当前块周边的像素。里面的像素点已经被癌细胞包围了,再放大就没有效果了。

为了进一步加快速度,我们在一次调用中对所有被选择沿相同方向扩展的像素执行膨胀。因此,每次迭代,我们最多执行 8 次扩张操作。

这使得代码相对较快(我测试了高达 1000x1000 的网格)。此外,它在所有迭代中保持相同的时间(不会随着网格开始填满而减慢)。

这是我的实现:

%# initial grid
img = false(500,500);

%# pick 10 random cells, and set them as cancerous
img(randi(numel(img),[10 1])) = true;

%# show initial image
hImg = imshow(img, 'Border','tight', 'InitialMag',100);

%# build all possible structing elements
%# each one dilates in one of the 8 possible directions
SE = repmat([0 0 0; 0 1 0; 0 0 0],[1 1 8]);
SE([1:4 6:9] + 9*(0:7)) = 1;

%# run simulation until all cells have cancer
BW = false(size(img));
while ~all(img(:)) && ishandle(hImg)
    %# find pixels on the perimeter of all "blocks"
    on = find(bwperim(img,8));

    %# percentage chance of division
    on = on( rand(size(on)) > 0.5 );    %# 50% probability of cell division
    if isempty(on), continue; end

    %# decide on a direction for each pixel
    d = randi(size(SE,3),[numel(on) 1]);

    %# group pixels according to direction chosen
    dd = accumarray(d, on, [8 1], @(x){x});

    %# dilate each group of pixels in the chosen directions
    %# to speed up, we perform one dilation for all pixels with same direction
    for i=1:8
        %# start with an image with only those pixels set
        BW(:) = false;
        BW(dd{i}) = true;

        %# dilate in the specified direction
        BW = imdilate(BW, SE(:,:,i));

        %# add results to final image
        img = img | BW;
    end

    %# show new image
    set(hImg, 'CData',img)
    drawnow
end

我还在 500x500 网格上创建了​​一个模拟动画,其中包含 10 个随机初始癌细胞(警告:.gif 图像大小约为 1MB,因此根据您的连接可能需要一些时间来加载)

模拟动画

于 2012-06-23T01:09:44.917 回答