scanpy.pl.tracksplot
- scanpy.pl.tracksplot(adata, var_names, groupby, use_raw=None, log=False, dendrogram=False, gene_symbols=None, var_group_positions=None, var_group_labels=None, layer=None, show=None, save=None, figsize=None, **kwds)
In this type of plot each var_name is plotted as a filled line plot where the y values correspond to the var_name values and x is each of the cells. Best results are obtained when using raw counts that are not log.
groupby
is required to sort and order the values using the respective group and should be a categorical value.- Parameters:
- adata :
AnnData
Annotated data matrix.
- var_names :
Union
[str
,Sequence
[str
],Mapping
[str
,Union
[str
,Sequence
[str
]]]] var_names
should be a valid subset ofadata.var_names
. Ifvar_names
is a mapping, then the key is used as label to group the values (seevar_group_labels
). The mapping values should be sequences of validadata.var_names
. In this case either coloring or ‘brackets’ are used for the grouping of var names depending on the plot. Whenvar_names
is a mapping, then thevar_group_labels
andvar_group_positions
are set.- groupby :
Union
[str
,Sequence
[str
]] The key of the observation grouping to consider.
- use_raw :
Optional
[bool
] (default:None
) Use
raw
attribute ofadata
if present.- log :
bool
(default:False
) Plot on logarithmic axis.
- num_categories
Only used if groupby observation is not categorical. This value determines the number of groups into which the groupby observation should be subdivided.
- categories_order
Order in which to show the categories. Note: add_dendrogram or add_totals can change the categories order.
- figsize :
Optional
[Tuple
[float
,float
]] (default:None
) Figure size when
multi_panel=True
. Otherwise thercParam['figure.figsize]
value is used. Format is (width, height)- dendrogram :
Union
[bool
,str
] (default:False
) If True or a valid dendrogram key, a dendrogram based on the hierarchical clustering between the
groupby
categories is added. The dendrogram information is computed usingscanpy.tl.dendrogram()
. Iftl.dendrogram
has not been called previously the function is called with default parameters.- gene_symbols :
Optional
[str
] (default:None
) Column name in
.var
DataFrame that stores gene symbols. By defaultvar_names
refer to the index column of the.var
DataFrame. Setting this option allows alternative names to be used.- var_group_positions :
Optional
[Sequence
[Tuple
[int
,int
]]] (default:None
) Use this parameter to highlight groups of
var_names
. This will draw a ‘bracket’ or a color block between the given start and end positions. If the parametervar_group_labels
is set, the corresponding labels are added on top/left. E.g.var_group_positions=[(4,10)]
will add a bracket between the fourthvar_name
and the tenthvar_name
. By giving more positions, more brackets/color blocks are drawn.- var_group_labels :
Optional
[Sequence
[str
]] (default:None
) Labels for each of the
var_group_positions
that want to be highlighted.- var_group_rotation
Label rotation degrees. By default, labels larger than 4 characters are rotated 90 degrees.
- layer :
Optional
[str
] (default:None
) Name of the AnnData object layer that wants to be plotted. By default adata.raw.X is plotted. If
use_raw=False
is set, thenadata.X
is plotted. Iflayer
is set to a valid layer name, then the layer is plotted.layer
takes precedence overuse_raw
.- show :
Optional
[bool
] (default:None
) Show the plot, do not return axis.
- save :
Union
[str
,bool
,None
] (default:None
) If
True
or astr
, save the figure. A string is appended to the default filename. Infer the filetype if ending on {'.pdf'
,'.png'
,'.svg'
}.- ax
A matplotlib axes object. Only works if plotting a single component.
- **kwds
Are passed to
heatmap()
.
- adata :
- Returns:
: A list of
Axes
.
Examples
>>> import scanpy as sc >>> adata = sc.datasets.pbmc68k_reduced() >>> markers = ['C1QA', 'PSAP', 'CD79A', 'CD79B', 'CST3', 'LYZ'] >>> sc.pl.tracksplot(adata, markers, 'bulk_labels', dendrogram=True)
Using var_names as dict:
>>> markers = {'T-cell': 'CD3D', 'B-cell': 'CD79A', 'myeloid': 'CST3'} >>> sc.pl.heatmap(adata, markers, groupby='bulk_labels', dendrogram=True)
See also
pl.rank_genes_groups_tracksplot
to plot marker genes identified using the
rank_genes_groups()
function.