## ABSTRACT

Analysis of natural scene statistics has been a powerful approach for understanding neural coding in the auditory and visual systems. In the field of somatosensation, it has been more challenging to quantify the natural tactile scene, in part because somatosensory signals are so tightly linked to the animal's movements. The present work takes a step towards quantifying the natural tactile scene for the rat vibrissal system by simulating rat whisking motions to systematically investigate the probabilities of whisker–object contact in naturalistic environments. The simulations permit an exhaustive search through the complete space of possible contact patterns, thereby allowing for the characterization of the patterns that would most likely occur during long sequences of natural exploratory behavior. We specifically quantified the probabilities of ‘concomitant contact’, that is, given that a particular whisker makes contact with a surface during a whisk, what is the probability that each of the other whiskers will also make contact with the surface during that whisk? Probabilities of concomitant contact were quantified in simulations that assumed increasingly naturalistic conditions: first, the space of all possible head poses; second, the space of behaviorally preferred head poses as measured experimentally; and third, common head poses in environments such as cages and burrows. As environments became more naturalistic, the probability distributions shifted from exhibiting a ‘row-wise’ structure to a more diagonal structure. Results also reveal that the rat appears to use motor strategies (e.g. head pitches) that generate contact patterns that are particularly well suited to extract information in the presence of uncertainty.

## INTRODUCTION

The sensory and motor systems of an animal species co-evolve with its nervous system, and all three evolve within the context of that species' particular ethological niche. In the fields of vision and audition, the receptive fields of central neurons are tuned to the statistics of the ‘natural scene’, that is, the stimuli that the animal is likely to encounter in its natural environment. Thus, neurons of primary visual cortex are tuned for the horizontal orientations and low spatial frequencies that dominate natural scenes (Vinje and Gallant, 2000; David et al., 2004), whereas neurons of higher-level auditory centers select for the sounds of the animal's vocalizations in terms of auditory objects rather than invariant acoustic features (Margoliash, 1983, 1986; Rauschecker et al., 1995; Ohlemiller et al., 1996; Nelken, 2004).

Analysis of natural scene statistics has proven to be a powerful approach for understanding neural coding within audition and vision, but has yet to be generally applied to the somatosensory system, although some work has been done on texture (Ritt et al., 2008; Manfredi et al., 2014). In the context of somatosensation, the way in which the animal's tactile sensors sample the environment will obviously play an important role in determining the tactile sensory data it will acquire. The statistics of the tactile scene are thus shaped by at least four species-specific attributes, all mutually-tuned through evolution: (1) the morphology and kinematics of the animal's sensors; (2) the receptive fields of the animal's neurons; (3) the motor strategies the animal adopts; (4) the stimuli the animal is likely to encounter within its particular ethological niche.

In the present work, we use the rat vibrissal system to examine these four tightly intertwined attributes. First, we quantify how whisking kinematics and the geometry of the array underlie the vibrissal contact patterns on a surface. Second, we investigate the extent to which the statistics of these contact patterns might mirror the receptive fields of central neurons in the vibrissal system. Third, we examine how the animal's behavioral strategies exploit the contact patterns. And finally, we describe how the contact patterns are matched (or mismatched) to the more complex environments of a rectangular cage and a typical rodent burrow.

## RESULTS

### Head pitch alters the trajectory of the whiskers in world-centered coordinates

We began by examining how whisking kinematics couple with array morphology to direct the whiskers through the environment. The overall motion of the whisker array was simulated during a typical protraction (see Materials and methods) using kinematic trajectories obtained from behaving animals (Knutsen et al., 2008). Although the primary direction of whisking motion is rostrocaudal, whisking trajectories also contain small components of elevation in the dorsoventral direction (Bermejo et al., 2002; Knutsen et al., 2008); the whisker also exhibits roll about its own axis (Knutsen et al., 2008).

Previous studies have shown that the elevation component of whisking causes the trajectories of the vibrissal tips to align more closely with the pitch of the head (Huet and Hartmann, 2014), but the trajectories of the whiskers in world coordinates have not yet been examined. We quantified the average direction in which the array points at each time step during a protraction. To obtain this measure, the direction vector of each vibrissa (from base to tip) was normalized to a length of one. These normalized vectors were averaged to obtain a single 3D vector that represents the average orientation and position of the array at every time step.

The average direction vector for the right vibrissal array over the course of a 60 deg protraction is shown in Fig. 1A. The figure shows that in head-centered coordinates, the average movement of the array is almost entirely ventrally directed up to ∼45 deg protraction, with a small dorsally directed motion during additional protraction. To confirm the robustness of these results, simulations were run with the parameter variations described in the Materials and methods. Although the trajectories show some variability, the overall motion over the protraction remains robust (Fig. 1B).

The protraction trajectory shown in Fig. 1A,B is invariant to head pitch because it is defined in head-centered coordinates. However, the effect of this trajectory in world-centered coordinates depends strongly on head pitch. Fig. 1C illustrates an example of the transformation between head-centered and world-centered coordinates: the head coordinates (black vectors) remain unchanged regardless of head pitch; in contrast, world coordinates (dark-red and blue vectors) change dramatically.

Four examples comparing trajectories in head-centered and world-centered coordinates are shown in the subplots of Fig. 1D, for head pitches of 60 deg, 0 deg, −25 deg and −70 deg; the left column illustrates the trajectories of the whisker tips in head-centered coordinates. Each figure is simply a rotation of the others. In contrast, the right column shows that the *z*-component of the average direction vector in world coordinates varies strongly with head pitch.

This effect of head pitch on trajectories in world coordinates is generalized in Fig. 1E. Head pitches above −10 deg (purple) always orient the array upwards, whereas head pitches below −44 deg (cyan) always orient the array toward the ground. Pitches between these values (shown in black) are initially oriented toward the ground, but begin orienting upward as protraction increases. These intermediate pitches correspond to those commonly observed during locomotion.

The importance of these coordinate transformations is twofold. First, the trajectories in world coordinates reveal that head pitch will significantly alter the direction of the reaction forces generated when the whiskers make contact with a surface. Second, the kinematic complexities make it much less evident which whiskers are most likely to make contact with any given surface. We explore this issue below.

### Array morphology and whisking kinematics determine which whiskers are most likely to contact a surface

The types of simulations shown in Fig. 1 allow us to establish the probability distributions that characterize whisker–object contact patterns during natural exploratory behavior. The space of all possible surfaces is large, so we initially consider surfaces that are relatively flat on the scale of the vibrissal array (radius of curvature ≥∼11 inches).

The premise of these simulations is that we imagine that over a rat's lifetime, it will encounter surfaces at all distances and orientations relative to its head. To simulate this ‘uniform’ prior distribution, we systematically varied the distance of the rat's head to the surface and the orientation of the rat's head relative to the surface. We used 61 different values of distance, 37 values of yaw and 37 values of pitch to obtain 83,509 different (distance, yaw, pitch) configurations, as shown in Fig. 2A,B and described in the Materials and methods. Head roll was omitted because it has no effect on whisker–surface contact statistics when averaged over the configuration space.

A whisk was simulated by protracting all vibrissae 60 deg from rest. A vibrissa was determined to have contacted the surface if it intersected any part of the surface. If, for a particular configuration, a whisker was already in contact with the surface at the start of the whisk, it was termed a ‘resting contact’. Otherwise, if a whisker came into contact with a surface at some point during its trajectory, it was termed a ‘whisking contact’. Whisking contacts can be further characterized by θ_{impact}, the angle between the rat's rostrocaudal midline and the line tangent to the whisker base when the whisker first impacts the object (Fig. 2C).

Unsurprisingly, these simulations demonstrated that longer vibrissae can make contact with surfaces in a larger number of configurations, simply because they have a larger reach. The longest, caudal-most vibrissae are in resting or whisking contact with surfaces in the largest number of configurations (Fig. 2D) and the relationship between arc length and the number of configurations that generate vibrissal–object contact is strong (*r*=0.929; Fig. 2E).

One dominant feature of rat exploratory behavior, however, is the tendency for the rat to orient towards objects with right–left symmetry (Milani et al., 1989; Benedetti, 1995; Hemelt and Keller, 2007; Cohen et al., 2008; Grant et al., 2012). Fig. 2F illustrates the number of whisking contacts for each vibrissa when the rat faces a flat surface symmetrically. Under these conditions, vibrissae near the array center (especially B2, B3, C2 and C3) are most likely to be in whisking contact. Thus, when only whisking contacts are considered, the relationship between arc length and the number of contacts (for yaw=0 deg) is much weaker (*r*=0.585; Fig. 2G).

This shift in contacts from the caudal to central vibrissae is an effect of both the right–left symmetry as well as the restriction to whisking contacts. The caudal vibrissae that dominate the overall number of contacts are often in resting contact when the head is turned to the side relative to the surface. Vibrissae near the center (not the front) of the array are the most likely to make whisking contact once the rat has oriented toward it.

The results of Fig. 2 illustrate the probabilities that individual whiskers make contact with a surface, but neurons at central levels of the vibrissal trigeminal system are generally characterized by multi-whisker receptive fields. In the next section, we therefore begin to quantify the joint probabilities of contact.

### A row-wise structure emerges from the probabilities of contact with flat surfaces

To start the characterization of joint probabilities of contact, we ask: given that a particular whisker makes contact with a surface at some point during a whisk, what are the probabilities that each of the other whiskers will also make contact with the surface during that whisk? We termed these joint contacts ‘concomitant contacts’. Note that the analysis of concomitant contacts does not imply anything about the timing of contact or the protraction angles at which contacts occur. When considered over both resting and whisking contacts, concomitant contact tends to be associated with nearest-neighbor vibrissae (Fig. 3A). An exception to this rule is seen for the rostral-most vibrissae. These vibrissae are so short that if they make contact with the surface, every other vibrissa is almost guaranteed to be also in contact.

Fig. 3B focuses on whisking-contacts: again, both row and column neighbors are more likely to be concomitant contact with the target vibrissa. However, the neighborhood of concomitantly contacting vibrissae is much more confined around the target. For about half of the whiskers, the caudal neighbor within the row exhibits the greatest likelihood of concomitant contact. Note that the overall color of the graphs shifts dramatically from purple (high probability) in Fig. 3A to green (low probability) in Fig. 3B. Although the total number of contacts is expected to decrease, it is not at all obvious that the probabilities should drop so dramatically when the analysis is limited to whisking contacts. The reason for this overall drop is that many configurations with resting contacts involve a large number of whiskers (e.g. 25/31 for high values of yaw). Because so many whiskers are in concomitant contact, a nearest-neighbor structure does not emerge from these configurations. In contrast, the configurations dominated by whisking contacts tend to involve far fewer whiskers (e.g. 5/31 for values of yaw near 0 deg). These configurations exhibit a clear nearest-neighbor structure. When combined, the resting contact configurations tend to dominate, reducing the nearest-neighbor structure in the contact pattern.

To examine the structure of concomitant contacts independent of whisker identity, a single probability distribution around the target vibrissa was created by averaging over all target whiskers. The average, shown in Fig. 3C, reveals a clear nearest-neighbor effect: vibrissae closest to the target show the greatest likelihood of concomitant contact. Interestingly, the row-neighbors show the greatest likelihood of concomitant contact. When a similar analysis is performed for the data in Fig. 3B (whisking contact only), the row-wise effect is even more pronounced (Fig. 3D). We emphasize that the averaged probability plots (Fig. 3C,D) are not proxies for every vibrissae and must be interpreted in light of individual plots (Fig. 3A,B). Rostral whiskers generally show a flatter probability structure than the identity-independent average, whereas the caudal vibrissae generally exhibit a more isolated and column-wise structure. These variations are consistent with previous results demonstrating that contact patterns are different for vibrissae in different regions of the array (Hobbs et al., 2015).

### Simulations predict the patterns of concomitant contact observed behaviorally

The analyses of Figs 2 and 3 quantified the contact probabilities within the assumption that the rat is equally likely to encounter a surface at each of the 83,509 possible configurations. During a particular behavior, however, the distribution of head poses relative to a surface will incorporate only a subset of poses. We therefore quantified whisker contacts as rats explored a flat, vertical glass sheet to identify the contact patterns associated with natural exploration patterns.

Naive rats perched on a T-maze and stretched across a gap to explore a vertical wall (see Materials and methods). The rat's head pose and the rostral-most and caudal-most whisker angles were tracked. We then simulated the vibrissal–object contact patterns that would be generated as the rat moved its head and whiskers through the same trajectory. Finally, we compared the simulated contact patterns with the contact patterns observed experimentally.

Qualitatively, the probabilities of concomitant contact predicted from simulation (Fig. 4A) are a good match for the behaviorally observed probabilities (Fig. 4C). Both simulated and observed behavior exhibit nearest-neighbor groupings, with the highest probability of concomitant contact occurring for nearest neighbors within the row. Intriguingly, both simulated and observed behavior suggest a diagonal structure in their identity-independent (average probability) plots (Fig. 4B,D). The most obvious difference between behavior and simulation is that the probability of concomitant contact drops off more sharply in the behavioral data than in simulation. This difference is observed as an increased number of blue and green (low-probability) contacts in Fig. 4C contrasted with larger regions of purple contacts (high probability) in Fig. 4A. Possible sources for this discrepancy are described in the Discussion.

Quantitatively, the match between simulation and behavior is remarkably good. Fig. 4E overlays the probabilities of concomitant contact for each vibrissa, for both simulation and behavior. The median absolute difference between behavior and simulation is shown above the target vibrissa's subplot. Across all vibrissa, the median absolute difference is 7.13±8.04% (s.d.).

The strong match between simulation and experimental data now enables two novel analyses, described in the next two sections. First, simulations can begin to predict the rat's behavioral strategies based on optimization criteria. Second, we can start to form predictions for the vibrissal contact patterns that will be generated in more complex environments, such as a cage or a burrow.

### The pitch of the rat's head can be predicted based on a strategy of maximizing the number of vibrissae in contact in the presence of uncertainty

To predict the rat's behavioral strategies based on optimization criteria, we first identified (in simulation) the head poses for which the number of whisking contacts with a surface was large. The set of poses in which 20 or more vibrissae (including both right and left sides) made whisking contact with the surface is shown in supplementary material Movie 1 and Fig. 5A.

Inspection of Fig. 5A shows that the number of whisking contacts is globally maximized for a distance of 4 mm, a yaw of 0 deg, and a negative pitch value, between 0 deg and −10 deg. It is also clear, however, that the contact distribution exhibits an angled structure, sloping upwards in pitch over distance. When summed over all values of distance, the largest number of whisking contacts is actually found near positive (rather than negative) pitch values. This value can be thought of as the expectation value of the number of contacts with a surface at an unknown distance, and is maximized near pitch=20 deg (Fig. 5B).

Thus Fig. 5A,B show that the rat should pitch its head near 20 deg to maximize whisking contact with a vertical wall, if distance is unknown. By extension, the rat should pitch its head to −70 deg to maximize contacts with the ground, at a 90 deg angle to the wall. This bimodal prediction from the simulation is shown in Fig. 5C, collapsed across yaws and distances.

We compared these predictions with the head pitches observed as naive rats stretched across a gap to explore a vertical wall. Fig. 5D shows the length of time, across all rats and all trials, in which the rat's head was pitched at each angle. The behavioral data shows peaks near 25 deg and −60 deg, corresponding closely to the predicted values of 20 deg and −70 deg. These are the head poses ideal for optimizing the number of whisking contacts with flat surfaces in the environment when yaw and distance are unknown. These results suggest that the rat pitches its head not to maximize contacts globally, but rather to maximize the expected number of contacts in the presence of uncertainty.

### Simulation of exploratory patterns in a cage

Finally, we can begin to construct predictions for the vibrissal contact patterns that will be generated in morecomplex environments. The first five numbered panels in Fig. 6A show frames from supplementary material Movie 2, which simulates contact patterns as the rat approaches and explores a wall. As the rat approaches the wall with its head pitched near −30 deg (Fig. 6A, panel 1), its ventral vibrissae remain in contact with the ground continuously. When it explores the ground–wall intersection (panel 2), an increased number of central-column vibrissae make contact. Note that the rostral-most vibrissae are not in contact with either surface, even though the central vibrissae make extensive contact. In panels 3–5 the rat rears its head from low to high and back, bringing many whiskers into contact both before and during whisking. With a pitch of −15 deg (panel 3), whisking contacts with the wall occur at large protraction angles (the E1 vibrissa contacts the ground). Surprisingly, none of the rostral-most or caudal-most vibrissae touch the wall, although the distance is only 5 mm. As the head pitches upwards (panel 4), the central and ventral vibrissae contact the wall at much lower θ_{impact} even during rest, reminiscent of contacts made with the ground. When the head returns to a pitch of 0 deg (panel 5), the rostral vibrissae make contact at very large values of θ_{impact}.

Fig. 6B,C generalize this typical exploration pattern by placing the wall and ground at different distances from the rat's head, (wall, 0–40 mm; ground, 15–40 mm; both in 5 mm increments) and by rotating the head through pitches ranging from −40 deg to 20 deg and yaws between −5 and 5 deg (both in 5 deg increments). For all contacts, the identity-independent average probability map is highly intuitive: the six vibrissae closest to the target show a high probability of concomitant contact, with a slight preference for the row-neighbors. As shown previously in Fig. 3, limiting the analysis to whisking-only produces neighborhoods of concomitant contact that are more concentrated around the target vibrissa.

Fig. 6D (panels 6–8) continues the sequence in supplementary material Movie 2, as the rat begins to turn to orient alongside the wall. The number of resting contacts increases dramatically and the contact patterns on the two sides become asymmetric. Fig. 6E,F generalize this effect by allowing the head to vary in distance from the wall and ground (wall, 0–40 mm; ground, 15–40 mm; both in 5 mm increments) and in orientation (pitch from −40 deg to 20 deg and yaws between 85 and 95 deg, both in 5 deg increments). In this case, a column-wise structure begins to emerge.

The column-wise structure occurs because at high yaw the caudal vibrissae make contact in a large number of configurations, and these configurations are dominated by resting contacts. In distinct contrast to the previous patterns, the nearest column-neighbors therefore have almost the same probability of concomitant contact as the nearest row-neighbors. We reiterate that these averaged probabilities are not proxies for individual vibrissa. The rostral-most vibrissae would show relatively flat probability structures for all neighbors, while the caudal vibrissae would show a stronger column-wise structure.

### Simulation of exploratory patterns in a burrow

Fig. 7A shows frames from supplementary material Movie 3, depicting the vibrissal contact patterns with the walls of a burrow. In these frames, the rat is simulated to be slightly off-center within the burrow in order to break right–left symmetry. A small-to-average burrow size is approximately 60 mm in diameter (Pisano and Storer, 1948; Calhoun, 1963). In a burrow of this size (Fig. 7A, panel 1), most whiskers are in resting contact. As the burrow's diameter increases to 80 mm (panel 2), the number of resting contacts for the right (near) array decreases only slightly, while additional resting contacts on the left (far) array become whisking contacts. For even larger diameters (perhaps corresponding to a nest within a burrow system), the right array experiences an increasing number of whisking contacts, but a subset of nine vibrissae remains in resting contact. Even for the left array, a subset of six vibrissae remains in resting contact with the ground.

This simulation reveals several key trends in the contact patterns. First, for average-sized burrows (60–100 mm), a large number of vibrissae are in resting contact; in particular, the ventral vibrissae maintain resting contact regardless of burrow size. Therefore, almost all information the rat can gather about the features of the burrow will occur as whiskers bend and deform continuously against burrow surfaces. Second, although whisking brings some additional vibrissae into contact with the wall, these whisking contacts occur at relatively small protraction angles (blue/green colors). Large (orange and red) θ_{impact} values are absent from these patterns: the rat gains information with little to no whisking.

As in previous analyses, the results for the burrow are generalized in Fig. 7B,C. In simulation, the burrow diameter was varied between 40 and 100 mm in 10 mm increments. The head was offset from the ground between 15 and 35 mm, offset to the left/right by up to 5 mm, and pitched between −40 deg and 20 deg in 5 deg increments.

Similar to results for the side-wall, contact patterns in the burrow are dominated by resting contacts, but in this case, a diagonal structure emerges in the probabilities of concomitant contact. In the burrow, many dorsal and rostral vibrissae never contact the walls; these also fall along a diagonal. The identity-independent map further emphasizes this trend: nearest neighbors within the row have the highest probability of concomitant contact, followed by diagonal (caudal–dorsal, rostral–ventral) neighbors, and then nearest neighbors within the column.

## DISCUSSION

### Descriptions of the tactile scene: environmental and vibrissal priors and posteriors

Because the sense of touch is so closely associated with movement, it may initially seem as though the tactile natural scene is inextricable from the way that an animal moves its sensors. From a theoretical standpoint, however, it is essential to develop a statistical description of the tactile scene independent of the sensor that samples it (e.g. a hand or a whisker array) and independent of how it is sampled. This statistical description corresponds to the ‘environmental prior’ – it is the distribution of the geometry of the environment, unbiased by any sampling strategy. A critical open theoretical question is which variable(s) should be chosen to characterize the environmental prior. Answering this question will require significantly more investigation; however, our intuition is that the characterization must relate to the local curvature of object's surfaces (van der Horst and Kappers, 2008; Hartmann, 2009; Schroeder and Hartmann, 2012).

Alternatively, the statistics of the environment can be described in terms of variables related to a particular sensor of interest. The present work, for example, focuses on developing the vibrissal prior: the whiskers are simulated to uniformly sample the environment and generate an unbiased prior in terms of whisker variables. This study has specifically focused on two geometric variables of interest: the probabilities of concomitant contact and the horizontal angle of vibrissal–surface contact (θ_{impact}). The incorporation of sampling strategies and biases of the rat during natural behavior modifies the statistics of the prior to produce the ‘vibrissal posterior’. It is this description that provides insight into how different behavioral strategies may alter the stimuli available to the rat for tactual perception. The present work's characterization of both the prior and posterior distributions during exploration of a flat wall is the first step towards a full characterization of the tactile scene.

The simulation results here are only the beginning. A complete description of the tactile scene must incorporate not only sensor geometry but also sensor mechanics. For example, whiskers (or fingertips) can rest gently on a surface or push against a surface. Future work will extend the geometric approach of the present manuscript to quantify the vibrissal prior and posterior in terms of the mechanical signals at the base of each whisker.

### Sensing at multiple scales: the spatial structure of contact depends critically on head pose

A key finding of the present study is that the basic spatial structure of the vibrissal prior depends critically on head pose. This result leads to an understanding of rat exploratory behavior in which head pose is the primary factor that determines the spatial structure of the vibrissal–surface contact pattern, that is, which vibrissae can make contact with a surface during any given whisk. Individual whisker kinematics then drive the temporal structure of contact, shaping when whiskers make contact with a surface and the resulting contact sequence. Thus, both head and whisker movements will be critical for determining the full spatiotemporal receptive field structure associated with the tactile natural scene.

Multiple scales of sensor movement, reflected as an interaction between head and whiskers, are not unique to the rat whisker system. A common theme seen across all sensory modalities is that low-mass sensors (e.g. eyes, pinnae and fingers) are located on larger-mass body structures (e.g. head and hands). The low-mass sensors permit control laws that can be close to kinematic, while inertial forces play a more important role in the movements of the larger structures. This mechanical ‘layering’ permits movements and sampling at multiple spatiotemporal scales.

In the context of the vibrissal system, whiskers are low mass and whisking motions are fundamentally rhythmic (Moore et al., 2013). A kinematic control law for whisking thus seems very plausible (Quist et al., 2014). In contrast, head movements are approximately ten times slower than whisking (Towal and Hartmann, 2006; Grant et al., 2009) and are more likely to involve significant inertial effects. Fig. 1 shows that even large alterations of array geometry and whisking kinematics (up to 15% variation in basepoint locations, and variations in roll and elevation within their established ranges) will have surprisingly small effects on the overall array trajectory. This suggests that regardless of how the rat chooses to whisk, key features of the contact patterns will remain largely invariant for a given head pose.

On a neurophysiological level, these results suggest that the receptive fields of whisker-sensitive neurons must be interpreted in light of the vestibular and proprioceptive signals informing the rat about head pose relative to a surface. Specifically, Fig. 1 illustrates that the 3D angle at which a whisker will contact a surface will depend strongly on head pitch. Given that neurons of the trigeminal system are strongly directionally tuned, these angles will have a large effect on the tactile signals generated.

### Sensing at multiple scales: the temporal structure of contact depends strongly on the kinematics of each whisk

We do not expect the simulations to be able to predict the timings of whisker contact on the surface. The simulation can predict only the patterns that could be generated by a given head pose relative to the surface, not whether that particular pattern will occur. Contact timing and sequence depend strongly on how the rat chooses to whisk.

For example, one might ask how contact probabilities would change if the rat did not protract through 60 deg, but instead protracted only ∼5 deg after first whisker contact. This type of simulation implicitly makes an assertion about contact timing, because whiskers do not generally protract uniformly or synchronously. Caudal and rostral whiskers often move out of phase, caudal whiskers often protract through larger angles than rostral whiskers, whiskers move at different velocities and whiskers start protraction from different set-points. At any point during the natural whisk cycle, each whisker might thus be at a different angle relative to its biomechanical rest. If a simulation began with all whiskers at their biomechanical rest positions, and then forced whiskers to stop exactly 5 deg after first whisker contact, very unnatural patterns of contact could be generated. The simulations we used examine contact probabilities completely independent of timing information. In other words, they capture the statistics of contact within the ‘reachable space’ of the array (Huet and Hartmann, 2014; Hobbs et al., 2015).

### Match between simulations and behavior

Support for the importance of head pose in determining spatial characteristics of contact is clearly seen in Fig. 4. In this figure, head pose and the overall whisk cycle capture key spatial features of the vibrissal–surface contact patterns. Given only head pose and position of rostral and caudal-most whiskers, the simulations can quite accurately (7.13% median error) predict the behaviorally observed patterns of concomitant contact.

Although the mystacial pad deforms during whisking (Wineski, 1983; Knutsen, 2015), previous studies have suggested that the effect of basepoint motion on the spatial structure of the contact patterns will be small (Hobbs et al., 2015). Several other factors, however, may help explain discrepancies observed between simulation and behavior.

First, the simulations assumed that all rows protract uniformly, which causes the probabilities to be more similar between rows in simulation than in behavior. This is likely to be the main source of the column-wise trend in percentage error observed in Fig. 4E: error is generally observed to decrease from dorsal to ventral. Second, the simulated positions of vibrissae in columns 1–5 were found by linearly interpolating the angles of the rostral-most and caudal-most whiskers. Because caudal vibrissae often lag the rostral (Carvell and Simons, 1990; Towal and Hartmann, 2008; Grant et al., 2009), this linear interpolation will tend to overestimate the probability of concomitant contact, particularly for rostral vibrissae. Third, simulations omitted the effect of head roll, which may partially account for why simulations under-predict the number of caudal–dorsal whiskers that will make contact (Fig. 4).

### Contact probabilities and receptive field structure

The simulations in the present work permit a systematic, exhaustive search through the space of possible contact patterns, thereby allowing for an estimation of the patterns that could occur during long sequences of natural exploratory behavior. These patterns would be extremely difficult to monitor and quantify in purely behavioral experiments. Specifically, the present study explores the patterns of concomitant contact in increasingly naturalistic conditions: the space of all possible configurations (Fig. 3), the behaviorally preferred configurations while exploring a surface (Fig. 4), and common environments such as cages (Fig. 6) and burrows (Fig. 7).

The space of all possible configurations (Fig. 3) shows that given a ‘principal whisker’ in contact with an object, its caudal neighbor within a row will tend to exhibit the greatest probability of also being in contact. This result echoes that of Carvell and Simons (1990), who showed that adjacent vibrissae are more likely to make concomitant contact during a texture discrimination task. Limiting the analysis to whisking contacts produces probability distributions that are even more selective for the caudal neighbor within a row (Fig. 3B). This structure is consistent with the row-wise receptive fields found in central brain structures including the trigeminal nuclei, sensory thalamus, and barrel cortex (Jacquin et al., 1989). In future work we anticipate analyzing these probability distributions subject to sparse coding constraints (Olshausen and Field, 1996) to more specifically predict receptive field characteristics.

The probability distributions of Fig. 3, however, assume that the rat is equally likely to explore all possible head configurations. During actual behavior, the rat chooses to explore particular configurations more than others. When contact probabilities were calculated as a rat explored a vertical glass wall, the neighborhood of increased probability grows to include secondary and tertiary neighbors, and a diagonal structure begins to emerge (Fig. 4). The diagonal structure is strongly enhanced in simulations of contact patterns generated in burrows (Fig. 7). The unique diagonal contact structure found for burrows, which comprise a large portion of the rat's natural habitat, is a strong argument for enrichment ‘tubes’ in cages. Exposure to these naturalistic contact patterns may be particularly important for the formation of ‘normal’ cortical circuitry, especially during the critical period (Hirsch and Spinelli, 1970).

### Exploratory strategies appear to maximize number of whiskers in contact in the presence of uncertainty

The diagonal and/or row-wise contact patterns that occur during exploration are a direct result of vibrissal array morphology, which has been tuned (jointly with the nervous system) on an evolutionary time scale. At the level of ontology, each rat is likely to exhibit motor behaviors that exploit this neuromechanical tuning. For example, rats are well known to orient towards surfaces in order to place a large number of whiskers in contact. One possible optimization strategy the rat might be using is to globally maximize the total number of whiskers in contact with the surface (Mitchinson et al., 2007; Grant et al., 2012). Fig. 5A shows, however, that achieving that global maximum would require the rat to place and maintain its head in a very specific pose relative to the surface.

Alternatively, the rat might compromise between the number of vibrissal contacts and the ability to vary the extent of its reach. This enables a large, although not globally maximal, number of whiskers to contact the surface even when its exact location is unknown. Across all distances, simulation results show that a head pitch near −70 deg maximizes the number of whisking contacts with the ground, while a head pitch near 20 deg produces the greatest number of whisking contacts with a vertical wall. Behavioral results show that rats preferred to explore with their heads near these poses. In a novel environment, the rats thus appear to prefer to explore in orientations which maximize the expected value of contacts with a surface at an unknown distance.

In conclusion, we found that the interplay between head and whisker morphology, whisking kinematics, head motion and the environment establishes the patterns and probabilities of vibrissal–object contact most relevant to behavior. The patterns of concomitant contact described in the present work are only one of several elements required to classify the natural tactile scene; a full description will also include the forces and moments at the base of all whiskers. This initial description of contact statistics has already provided insight into the rat's behavior, as well as the probability structures to which neurons are likely to be tuned.

## MATERIALS AND METHODS

Many methods are adapted from our previous work (Hobbs et al., 2015).

### Modeling whisking kinematics in three dimensions

A three-dimensional (3D) model of the rat head and vibrissal array (Towal et al., 2011) was used to simulate contact patterns. Three variables quantify the orientation of the vibrissae: protraction angle (θ), elevation angle (φ) and the rotation about the vibrissa's own axis (ζ). Knutsen et al. (2008) provides slope values relating dθ and dφ, and dθ and dζ during natural whisking motion, which allows us to write the kinematic equations shown in Table 1. Inclusion of the three kinematic degrees of freedom (protraction, elevation and roll) captures the full three-dimensional trajectories of the whiskers (Bermejo et al., 2002; Knutsen et al., 2008).

### Simulation of surfaces of interest

When simulating the entire configuration space, the forward kinematics of the 31 vibrissae on the right mystacial pad were simulated by stepping in 0.1 deg increments through a maximum protraction of 60 deg from biomechanical rest (θ_{rest}). The entire configuration space was simulated as follows: the distance between the nose and the wall was increased between 0 and 60 mm in 1 mm increments; head pitch was varied between −90 deg and +90 deg in 5 deg increments; and head yaw was varied between −90 deg and +90 deg in 5 deg increments. Running the simulation over the full range of distance, yaw and pitch generates 83,509 (61×37×37) configurations.

A 60 deg protraction produces angle values (relative to the rostral-caudal midline of the body) between 98 deg and 164 deg for the various whiskers, consistent with the range of angles seen in other studies (Carvell and Simons, 1990; Gao et al., 2001; Hill et al., 2008). Data were then reflected about the midline to obtain the trajectories for the left vibrissal array. The burrow was simulated as an infinite cylinder centered on the *y*-axis. The cylinder's diameter was varied between 60 mm and 140 mm, in 20 mm increments. To break right–left symmetry, the rat's nose was simulated to be 5 mm to the right of the cylinder's central axis. The head was pitched downward at −40 deg and the nose was placed 20 mm above the burrow floor (Thé et al., 2013).

### Definitions of θ_{impact} and contact type

Each vibrissa was simulated with 100 nodes; a vibrissa was determined to have collided with a surface if the coordinates of any node intersected the surface boundary. The θ angle between the rat's rostrocaudal midline and the tangent to the vibrissal base at the instant of first vibrissa contact is defined as θ_{impact} (cf. Fig. 2C).

### Sensitivity analysis

Fig. 1B makes use of a sensitivity analysis described in depth in a previous study, showing that trends in θ_{impact} are robust to a wide variety of parameter changes (Hobbs et al., 2015). The sensitivity analysis varied the following parameters independently: whisker length was allowed to decrease by 20%; resting parameters φ_{rest} and ζ_{rest} were uniformly distributed between ±20% of nominal; kinematic relationships for φ and **ζ** were normally distributed within the error bounds listed in Table 1; basepoint location was varied up to 15% between neighboring basepoints and uniformly in any direction. A Monte Carlo analysis was performed with 42,000 different combinations of these altered parameter values. The average direction vector of the array was calculated for each of these 42,000; a randomly selected subset of 1000 trials is shown in Fig. 1B.

### Measurement of head orientation during exploration

All experimental work was approved in advance by Northwestern University's Animal Care and Use Committee. Simulation results were compared with results from a behavioral experiment performed under infrared lighting, in which five naive Long Evans rats (*Rattus* *norvegicus* Berkenhout 1769; female; 2.5 months) perched on a ledge and freely explored a vertical wall located on the other side of a gap. The wall was glass with no textural features and so large that the rats could not touch its edges. The first 1–3 s of behavior were recorded as the rats ‘gained a first impression’ of the novel surface, for a total of 86 whisks (8.61 s) across all rats.

Two high-speed video cameras (1000 frames per second, Photron FastCam 1024-PCI) were used to monitor the rats as they explored the wall. One camera captured a ‘bird's-eye’ view and the second a ‘head-on’ view through the glass wall. The rat's eyes and nose were tracked in both camera views, and the trajectories filtered at 20 Hz. The 3D head position and orientation were reconstructed from these 2D trajectories using standard image processing techniques (Hartley and Zisserman, 2003).

In simulation, a head pitch of 0 deg corresponds to the head angle that places the average row-planes of the whisker basepoints parallel to the ground (Towal et al., 2011). In the behavioral experiment, data were taken so that a 0 deg head pitch corresponds to the head angle at which the plane containing the nose and eyes is parallel to the *x*–*y* plane. The offset (about ∼19 deg) was added to the behavioral head pitches to match the definition used in simulation.

The front and back vibrissae on both sides of the face were tracked in the top-down camera view. The effect of head pitch on the measured (raw) protraction angles was corrected for by projecting the angles into the *x*–*y* plane. Protraction angles were filtered at 100 Hz to remove noise. The protraction angles for the interior columns of the array were linearly interpolated between the measured front and back angles.

Vibrissal–object contacts were detected and identified using methods described previously (Towal and Hartmann, 2010). A series of optics produced a 2 mm plane of collimated laser light (975 nm) directly in front of the glass sheet. When whiskers contacted the glass, they interrupted the planar beam, scattering light at the points of contact. High-speed video cameras captured these points of light. Any points of light below the intensity threshold (established by independent image calibration) were discarded. All other contact points were analyzed using image processing algorithms based on particle tracking velocimetry.

## Acknowledgements

We thank Sara A. Solla for many helpful discussions and comments.

## Footnotes

**Author contributions**

J.A.H. and R.B.T. performed the experiments. J.A.H., R.B.T. and M.J.Z.H. analyzed the data. J.A.H. and M.J.Z.H. prepared the manuscript and R.B.T. helped to edit the manuscript.

**Funding**

This work was supported by National Science Foundation awards [IOS-0818414, IOS-0846088 and EFRI-0938007 to M.J.Z.H.].

## References

**Competing interests**

The authors declare no competing or financial interests.