How do I complete my matlab code for a given formula.?

% This is done through k-Means clustering and binarization to derive a 
% binary classification of each parcel for cropland and non-cropland pixels.
% This process results in 2 clusters per parcel object for which the 
% mean values are calculated. Herein, a threshold for pixels higher 
% than 0.4 was defined for the derivation of CAF:
% 𝐶𝐴𝐹 = Σ𝑝 > 0 . 4 / Σ𝑝
% where Σp is the sum of all pixels and Σ𝑝 > 0 . 4 the sum of pixels
% which are located within the cluster with a mean NDVI higher than 0.4.
% open NDVI file
img = imread('ndvi.tif');
% pickup the shape files of parcels
d = uigetdir(pwd, 'Select a folder');
shapefiles = dir(fullfile(d, '*.shp'));
for m = 1:length(shapefiles)
    shapefile = shapefiles(m);
    disp(shapefile.name);
    S = shaperead(shapefile.name);
    polygon = polyshape([S.X], [S.Y]);
    % Create a logical mask
    logical_mask = inpolygon(lon, lat, polygon.Vertices(:, 1), polygon.Vertices(:, 2));
    % Use the logical mask to extract data from ndvi image for parcel
    parcel_ndvi = img(logical_mask);
    % Apply k-means clustering with k=2 
    k = 2; 
    [idx, centroids] = kmeans(parcel_ndvi, k); 
    X = mean(idx(:,1));
    Y = mean(idx(:,2));
    
end

I would appreciate any help to complete above code as per comments given at the top of this code.

Deve

[Steven Lord]

It's not clear to me what help you're looking for. Please ask a specific question about where you're having difficulty and we may be able to provide some guidance about that specific area in your code.

[Devendra]

How to find cropland and non cropland cluster (which is which out of these two clusters) using the following matlab code line

% Apply k-means clustering with k=2 
k = 2; 
[idx, centroids] = kmeans(parcel_ndvi, k); 

thanks.

Dave

[Walter Roberson]

You could take

mean(parcel_ndvi(idx == 1))

and compare it to

mean(parcel_ndvi(idx == 2))

to try to decide which cluster is which ?

[Devendra]

if mean(parcel_ndvi(idx == 1)) > mean(parcel_ndvi(idx == 2)) then

cropland = idx(1)

noncropland = idx(2)

otherwise

cropland = idx(2)

noncropland = idx(1)

endif

How do I write above lines in matlab code ?

Thanks a lot for your kind help.

Dave

[Walter Roberson]

if mean(parcel_ndvi(idx == 1)) > mean(parcel_ndvi(idx == 2))
    cropland = 1;
    noncropland = 2;
else
    cropland = 2;
    noncropland = 1;
end

Answer 1

[Walter Roberson]

The first thing you have to do is figure out which of the two clusters corresponds to which condition.

kmeans() clustering with two clusters can return either index being the one of interest -- even if you give it options for initial cluster centers, because of the way that cluster updates go, it is possible for the two clusters to effectively exchange identities.

Once you have figured out which is which, you would

mask = idx == CLUSTER_OF_INTEREST;
CAF = nnz(parcel_ndvi(mask) > 0.4) / (size(img,1)*size(img,2))

[Devendra]

mask = idx == CLUSTER_OF_INTEREST;

CAF = nnz(parcel_ndvi(mask) > 0.4) / (size(img,1)*size(img,2))

it should be (size(idx,1)*size(idx,2)) or (size(img,1)*size(img,2)) ??

Thanks for your kind help.

Dave

[Walter Roberson]

size(img,2)*size(img,2) is the number of pixels covered by img

[Devendra]

Thank you very much for your kind help. I am using following matlab lines to create mask over cropland.

if mean(parcel_ndvi(idx == 1)) > mean(parcel_ndvi(idx == 2))

cropland = 1;

noncropland = 2;

else

cropland = 2;

noncropland = 1;

end

mask = idx == cropland

I want to mask only those pixels under cropland whose ndvi values are greater than 0.4. I request you to please suggest me how to create the said mask using matlab code. I would appreciate your kind help.

Devendra