Skip to content

gate_finder

Short Description

sm.pl.gate_finder: The function opens the OME-TIFF image inside Napari and overlays points to help with the identifying manual gates for each marker. Use the sm.pp.rescale function to apply the identified gates to your data.

Function

gate_finder(image_path, adata, marker_of_interest, layer='raw', log=True, from_gate=6, to_gate=8, increment=0.1, markers=None, channel_names='default', flip_y=True, x_coordinate='X_centroid', y_coordinate='Y_centroid', point_size=10, imageid='imageid', subset=None, seg_mask=None, **kwargs)

Parameters:

Name Type Description Default
image_path str

Location to the image file.

 required
adata

Ann Data Object

required
marker_of_interest str

Marker for which gate is to be defined e.g. 'CD45'.

 required
layer str

The layer in adata.layers that contains the expression data to gate. If None, adata.X is used. use raw to use the data stored in adata.raw.X

 'raw'
log bool

Log transform the data before gating.

 True
from_gate int

Start value gate of interest.

 6
to_gate int

End value of the gate of interest.

 8
flip_y bool

Flip the Y-axis if needed. Some algorithms output the XY with the Y-coordinates flipped. If the image overlays do not align to the cells, try again by setting this to False.

 True
increment float

Increments between the start and end values.

 0.1
markers str

Additional markers to be included in the plot for evaluation.

 None
channel_names list

List of channels in the image in the exact order as image. The default is adata.uns['all_markers']

 'default'
x_coordinate str

X axis coordinate column name in AnnData object.

 'X_centroid'
y_coordinate str

Y axis coordinate column name in AnnData object.

 'Y_centroid'
point_size int

point size in the napari plot.

 10
imageid str

Column name of the column containing the image id.

 'imageid'
subset str

imageid of a single image to be subsetted for analyis.

 None
seg_mask str

Location to the segmentation mask file.

 None 
1
2
3
4
5
6
    image_path = '/Users/aj/Desktop/ptcl_tma/image.ome.tif'
    sm.pl.gate_finder (image_path, adata, marker_of_interest='CD45',
                 from_gate = 6, to_gate = 8, increment = 0.1,
                 markers=['DNA10'], channel_names = 'default',
                 x_coordinate='X_position',y_coordinate='Y_position',point_size=10,
                 subset= '77', seg_mask=None)
Source code in scimap/plotting/_gate_finder.py 
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
def gate_finder (image_path, adata, marker_of_interest, layer='raw', log=True,
                 from_gate = 6, to_gate = 8, increment = 0.1,
                 markers=None, channel_names = 'default', flip_y=True,
                 x_coordinate='X_centroid',y_coordinate='Y_centroid',
                 point_size=10,imageid='imageid',subset=None,seg_mask=None,**kwargs):
    """
Parameters:
    image_path (str):  
        Location to the image file.

    adata : Ann Data Object  

    marker_of_interest (str):  
        Marker for which gate is to be defined e.g. 'CD45'.

    layer (str): 
        The layer in adata.layers that contains the expression data to gate. 
        If None, adata.X is used. use `raw` to use the data stored in `adata.raw.X`

    log (bool):  
        Log transform the data before gating.

    from_gate (int):   
        Start value gate of interest.

    to_gate (int):    
        End value of the gate of interest.

    flip_y (bool):  
        Flip the Y-axis if needed. Some algorithms output the XY with the Y-coordinates flipped.
        If the image overlays do not align to the cells, try again by setting this to `False`.

    increment (float):  
        Increments between the start and end values.

    markers (str):  
        Additional markers to be included in the plot for evaluation.

    channel_names (list):  
        List of channels in the image in the exact order as image. The default is `adata.uns['all_markers']`

    x_coordinate (str):  
        X axis coordinate column name in AnnData object.

    y_coordinate (str):  
        Y axis coordinate column name in AnnData object.

    point_size (int):  
        point size in the napari plot.

    imageid (str):  
        Column name of the column containing the image id.

    subset (str):  
        imageid of a single image to be subsetted for analyis.

    seg_mask (str):  
        Location to the segmentation mask file.

    **kwargs  
        Other arguments that can be passed to napari viewer.

Example:
```python
    image_path = '/Users/aj/Desktop/ptcl_tma/image.ome.tif'
    sm.pl.gate_finder (image_path, adata, marker_of_interest='CD45',
                 from_gate = 6, to_gate = 8, increment = 0.1,
                 markers=['DNA10'], channel_names = 'default',
                 x_coordinate='X_position',y_coordinate='Y_position',point_size=10,
                 subset= '77', seg_mask=None)
```
    """

    # If no raw data is available make a copy
    if adata.raw is None:
        adata.raw = adata

    # subset data if neede
    if subset is not None:
        if isinstance (subset, str):
            subset = [subset]
        if layer == 'raw':
            bdata=adata.copy()
            bdata.X = adata.raw.X
            bdata = bdata[bdata.obs[imageid].isin(subset)]
        else:
            bdata=adata.copy()
            bdata = bdata[bdata.obs[imageid].isin(subset)]
    else:
        bdata=adata.copy()

    # isolate the data
    if layer is None:
        data = pd.DataFrame(bdata.X, index=bdata.obs.index, columns=bdata.var.index)[[marker_of_interest]]
    elif layer == 'raw':
        data = pd.DataFrame(bdata.raw.X, index=bdata.obs.index, columns=bdata.var.index)[[marker_of_interest]]
    else:
        data = pd.DataFrame(bdata.layers[layer], index=bdata.obs.index, columns=bdata.var.index)[[marker_of_interest]]

    if log is True:
        data = np.log1p(data)


    # Copy of the raw data if it exisits
    #if adata.raw is not None:
    #    adata.X = adata.raw.X

    # Plot only the Image that is requested
    #if subset is not None:
    #    adata = adata[adata.obs[imageid] == subset]

    # Make a copy of the data with the marker of interest
    #data = pd.DataFrame(np.log1p(adata.X), columns = adata.var.index, index= adata.obs.index)[[marker_of_interest]]

    # Generate a dataframe with various gates
    def gate (g, d):
        dd = d.values
        dd = np.where(dd < g, np.nan, dd)
        np.warnings.filterwarnings('ignore')
        dd = np.where(dd > g, 1, dd)
        dd = pd.DataFrame(dd, index = d.index, columns = ['gate-' + str(g)])
        return dd

    # Identify the list of increments
    inc = list(np.arange (from_gate, to_gate, increment))
    inc = [round(num,3) for num in inc]

    # Apply the function
    r_gate = lambda x: gate(g=x, d=data) # Create lamda function
    gated_data = list(map(r_gate, inc)) # Apply function
    # Concat all the results into a single dataframe
    gates = pd.concat(gated_data, axis=1)


    # Recover the channel names from adata
    if channel_names == 'default':
        channel_names = adata.uns['all_markers']
    else:
        channel_names = channel_names

    # if markers is a string convert to list
    if isinstance(markers, str):
        markers = [markers]

    # Index of the marker of interest and corresponding names
    if markers is not None:
        markers.extend([marker_of_interest])
        idx = np.where(np.isin(channel_names,markers))[0]
        channel_names = [channel_names[i] for i in idx]
    else:
        idx = list(range(len(channel_names)))
        channel_names = channel_names


    # Load the segmentation mask
    if seg_mask is not None:
        seg_m = tiff.imread(seg_mask)
        if (len(seg_m.shape) > 2) and (seg_m.shape[0] > 1):
            seg_m = seg_m[0]


    ##########################################################################
    # Visulaisation using Napari

    # load OME TIFF
    if os.path.isfile(image_path) is True: 
        # Load the image    
        image = tiff.TiffFile(image_path, is_ome=False)
        z = zarr.open(image.aszarr(), mode='r') # convert image to Zarr array
        # Identify the number of pyramids and number of channels
        n_levels = len(image.series[0].levels) # pyramid  
        # If and if not pyramids are available
        if n_levels > 1:
            pyramid = [da.from_zarr(z[i]) for i in range(n_levels)]
            multiscale = True
        else:
            pyramid = da.from_zarr(z)
            multiscale = False   
        # subset channels of interest
        if markers is not None:
            if n_levels > 1:
                for i in range(n_levels-1):
                    pyramid[i] = pyramid[i][idx, :, :]
                n_channels = pyramid[0].shape[0] # identify the number of channels
            else:
                pyramid = pyramid[idx, :, :]
                n_channels = pyramid.shape[0] # identify the number of channels
        else:
            if n_levels > 1:
                n_channels = pyramid[0].shape[0]
            else:
                n_channels = pyramid.shape[0]\
        # check if channel names have been passed to all channels
        if channel_names is not None:
            assert n_channels == len(channel_names), (
                f'number of channel names ({len(channel_names)}) must '
                f'match number of channels ({n_channels})'
            )

        # Load the viewer
        viewer = napari.view_image(
        pyramid,
        channel_axis = 0,
        multiscale=multiscale,
        name = None if channel_names is None else channel_names,
        visible = False, **kwargs)


    # Operations on the ZARR image
    # check the format of image
    if os.path.isfile(image_path) is False: 
        #print(image_path)
        viewer = napari.Viewer()
        viewer.open(image_path, multiscale=True,
                    visible=False,
                    name = None if channel_names is None else channel_names)

    # Add the seg mask
    if seg_mask is not None:
        viewer.add_labels(seg_m, name='segmentation mask', visible=False)

    # subset the gates to include only the image of interest
    gates = gates.loc[bdata.obs.index,]

    # Add gating layer
    def add_phenotype_layer (adata, gates, phenotype_layer,x,y,viewer,point_size):
        cells = gates[gates[phenotype_layer] == 1].index
        coordinates = adata[cells]
        # Flip Y axis if needed
        if flip_y is True:
            coordinates = pd.DataFrame({'y': coordinates.obs[y],'x': coordinates.obs[x]})
        else:  
            coordinates = pd.DataFrame({'x': coordinates.obs[x],'y': coordinates.obs[y]})
        #points = coordinates.values.tolist()
        points = coordinates.values
        #import time
        #start = time.time()
        viewer.add_points(points, size=point_size, face_color='white',visible=False,name=phenotype_layer)
        #stop = time.time()
        #print(stop-start)


    # Run the function on all gating layer
    for i in gates.columns:
        add_phenotype_layer (adata=bdata, gates=gates, 
                             phenotype_layer=i, x=x_coordinate, y=y_coordinate, 
                             viewer=viewer, point_size=point_size)