curtains

Detecting DNA curtains

The smtools.curtains module is designed to locate, and fit the intensity profile of individual DNA molecules within a DNA curtain. For a description of how it works see our walkthroughs part 1, part 2, and part 3.

import smtools.curtains as cs
import smtools.alignment as al
from scipy.ndimage.interpolation import rotate
import matplotlib.pyplot as plt
import smtools.testdata as test
im = test.test_curtain()
ch1,ch2 = al.im_split(im)
angle = cs.find_rotang(ch2)
rotated_ch2 = rotate(ch2,angle)
bounds, mask = cs.find_curtain(rotated_ch2)
strands = cs.find_DNA(rotated_ch2,bounds)
DNAs = cs.fit_DNA(rotated_ch2, strands)
fig = plt.figure(figsize=(5,10))
plt.axis('off')
plt.imshow(rotated_ch2)
for x0,x1,y in DNAs:
    plt.plot([x0,x1],[y,y],"r.", markersize = 5)
plt.show()
_images/curtain_finder.png

The curtains module

curtains.find_rotang(Image, line_length=40, theta=None, line_gap=5, tilt=0.1)

Finds the general rotation of an image containing DNA curtains. Uses a Hough line Transform to approximate DNA locations. Returns the average angle of all hough lines.

Parameters:
  • Image – 2D image array
  • line_length

    int, minimum accepted length of line detected by the Hough line Transform Default is 40.

  • theta

    array, angles at which to compute the Hough line Transform in radians. Default is None

  • line_gap – int, maximum gap between pixels allowed when forming a line. default is 5
  • tilt

    int, when theta is None, the tilt value defines the range of thetas to be used in the Hough line Transform as (pi/2 - tilt:pi/2 + tilt). Defaults to 0.1 radians

Returns:

measure rotation angle in radians:

Example:
>>> import smtools.curtain as cs
>>> import smtools.testdata as test
>>> im = test.test_curtain()
>>> ch1,ch2 = al.im_split(i)
>>> print(cs.find_rotang(ch2))
0.6936316728744394
curtains.find_curtain(Image, distance=50, line_length=40, theta=None, line_gap=5, tilt=0.1, window=15, order=3, maxline=70)
Parameters:
  • Image – 2D image array
  • distance – int, passes to peak_finder. defines minimum spacing between two separate curtains. Default is 50.
  • line_length

    int, minimum accepted length of line detected by the Hough line Transform Default is 40.

  • theta

    array, angles at which to compute the Hough line Transform in radians. Default is None

  • line_gap – int, maximum gap between pixels allowed when forming a line. default is 5
  • tilt

    int, when theta is None, the tilt value defines the range of thetas to be used in the Hough line Transform as (pi/2 - tilt:pi/2 + tilt). Defaults to 0.1 radians

  • window – int, window value passed to the savgol_filter
  • order

    int, order value passed to the savgol_filter

  • maxline

    int, maximum size line detected by the Hough line Transform. defaults to 70.

Returns:

List of tuples containing (x_min, x_max, y_min, y_max) for each curtain in the image.

Returns:

Boolean mask of the image. Evaluates True within the bounds of the detected curtains, and False outside.

Example:
>>> import smtools.curtain as cs
>>> import smtools.testdata as test
>>> from scipy.ndimage.interpolation import rotate
>>> im = test.test_curtain()
>>> ch1,ch2 = al.im_split(i)
>>> angle = cs.find_rotang(ch2)
>>> rotated_ch2 = rotate(ch2,angle)
>>> bounds, mask = cs.find_curtain(rotated_ch2)
>>> print(bounds)
[(1, 76, 0, 515), (89, 157, 0, 513), (167, 235, 0, 515)]
curtains.find_DNA(Image, Bounds, prominence=None)
Parameters:
  • Image – 2D image array
  • Bounds – 1D List of tuples containing (x_min, x_max, y_min, y_max) for each curtain in the image. written to accept output from toolbox.curtains.find_curtain
  • prominence

    passes to peak_finder Defaults to (20, 1.5 * max(array)

Returns:

List of tuples: (top, bottom, and center). Locations of each DNA strand detected in pixels.

Example:
>>> import smtools.curtain as cs
>>> import smtools.testdata as test
>>> from scipy.ndimage.interpolation import rotate
>>> im = test.test_curtain()
>>> ch1,ch2 = al.im_split(i)
>>> angle = cs.find_rotang(ch2)
>>> rotated_ch2 = rotate(ch2,angle)
>>> bounds, mask = cs.find_curtain(rotated_ch2)
>>> strands = cs.find_DNA(rotated_ch2,bounds)
>>> print(len(strands))
269
curtains.fit_DNA(Image, locations, constraints=None, init_guess_length=None, init_guess_center=None, pad=2, factor=1.3, min_length=45, cen_err=0.3, cen_dev=3.0)
Parameters:
  • Image – 2D image array
  • locations – 1List of tuples: (top, bottom, and center). Locations of each DNA strand detected in pixels. written to accept output from toolbox.curtains.find_DNA
  • constraints – 1D list, size 5. bounds on the call to curve_fit when fitting toolbox.misc.super_gauss_function. see bounds Default is ([0, 0, 0, 0, 5.],[np.inf, np.inf, np.inf, np.inf, np.inf]).
  • init_guess_length – 1D list, size 5. initial guess for the parameters for the call to curve_fit when fitting toolbox.misc.super_gauss_function. see p0 Default is [.2, 1., np.median(xdata), 20, 6].
  • init_guess_center

    1D list, size 5. initial guess for the parameters for the call to curve_fit when fitting toolbox.misc.mundane_gauss_function. see p0 Default is [1, 1, 2, 1].

  • pad – int, number of adjacent rows of pixels to grab for the fit, Default is 2 pixels out on each side.
  • factor – float, multiplier used in end detection. defines the intensity threshold at the ends as maximum - ((maximum - minimum)/factor)
  • min_length – int, minimum allowed length of DNA. Default is 45.
  • cen_err – float, maximum error ratio on the mean for the call to curve_fit when fitting toolbox.misc.mundane_gauss_function. Defined as np.sqrt(pcov[2][2])/popt[2] Defualts to 0.3.
  • cen_dev

    float, maximum allowed standard deviation for the call to curve_fit when fitting toolbox.misc.mundane_gauss_function. Defualts to 3.

Returns:

List of tuples: (top, bottom, and center). Locations of each DNA strand detected as defined by the fitting with toolbox.misc.super_gauss_function and toolbox.misc.mundane_gauss_function.

Example:
>>> import smtools.curtain as cs
>>> import smtools.testdata as test
>>> from scipy.ndimage.interpolation import rotate
>>> im = test.test_curtain()
>>> ch1,ch2 = al.im_split(i)
>>> angle = cs.find_rotang(ch2)
>>> rotated_ch2 = rotate(ch2,angle)
>>> bounds, mask = cs.find_curtain(rotated_ch2)
>>> strands = cs.find_DNA(rotated_ch2,bounds)
>>> DNAs = cs.fit_DNA(rotated_ch2, strands)
>>> print(len(DNAs))
200
>>> import matplotlib.pyplot as plt
>>> plt.figure()
>>> plt.imshow(rotated_ch2)
>>> for x0,x1,y in DNAs:
>>>     plt.plot([x0,x1],[y,y],"r.", markersize = 5)
>>> plt.show()
class curtains.Curtain(DNA)

Collection into a Curtain .

DNA : list
Accepts output from fit DNA
points : list or list of tuples
points to apply to DNA
>>> import smtools.curtain as cs
>>> import smtools.testdata as test
>>> from scipy.ndimage.interpolation import rotate
>>> im = test.test_curtain()
>>> ch1,ch2 = al.im_split(i)
>>> angle = cs.find_rotang(ch2)
>>> rotated_ch2 = rotate(ch2,angle)
>>> bounds, mask = cs.find_curtain(rotated_ch2)
>>> strands = cs.find_DNA(rotated_ch2,bounds)
>>> DNAs = cs.fit_DNA(rotated_ch2, strands)
>>> strands = Curtain(DNAs)
>>> strands.assign_points(test.test_points())
>>> data = strands.get_DNA()
>>> print(data)
add_point(key, point)
Parameters:
  • key – keys to DNA_dictionary
  • point – (x,y)
assign_points(points, frame=None, pad=5, max_offset=1)
Parameters:
  • points – 1D list of (x,y) points
  • frame – frame number if data are from image stack
  • pad – additional extension to length of DNA applied to top and bottom position when assigning points.
  • max_offset – maximum distance from center of DNA strand that still allows a point to be assigned to it.
assign_point_stack(points_list, pad=5, max_offset=1)
Parameters:
  • points_list – list of 1D lists of (x,y) points
  • pad – additional extension to length of DNA applied to top and bottom position when assigning points.
  • max_offset – max_offset: maximum distance from center of DNA strand that still allows a point to be assigned to it.
pairwise()
Returns:
per_strand_counts()
Returns:
get_DNA()
Returns:dict, Keys are (top, bottom, center) for each DNA strand. Entries are lists of points (x, y) or (x, y, frame)