Xplor-NIH (frequently) Asked Questions

1. How do I generate parameters for a small molecule?
   
2. How do I make a cis-residue (such as cis-proile)?
   
3. script to add-in atoms whose positions are initially unknown
   
4. using the IVM (internal variable modules) for torsion angle dynamics and more
   
5. I see errors like: Powell::step: irregular exit: Line search abandoned: gradient may be incorrect.
   
6. Radius of gyration potential term
   
7. DELPHIC torsion database of phi-psi angles
   
8. DELPHIC torsion database of phi-psi angles: potential definition.
   
9. How do I determin the structure of a homo-dimer?
   
10. dipolar couplings: the definition of the coefficients
    
11. dipolar couplings: how are force constants weighted between different experiments
    
12. Strange characters in output files or binary output files.
    
13. Xplor-NIH for similar computers.
    
14. What are random number seeds?
    
15. why do I get the error: %ASSFIL-ERR: no free unit available?
    
16. Why do I get the error %COPYST-ERR: ST2MAX too small. Check input file?
    
17. Why do I get the error cannot restore segment prot after reloc: Permission denied
    
18. Are Intel-based Macs supported?
    
19. Can Xplor-NIH 2.24 (and maybe earlier) run on Maxc OS X 10.6?
    
20. Why doesn't Xplor-NIH run on fedora core 5?
    
21. does Xplor-NIH support parallel computer architectures?
    
22. How do I run Xplor-NIH on a cluster using the PBS queuing system?
    

You might try the one of these resources:

1) Check the PDB's ligand Expo: http://ligand-expo.rcsb.org/.
   If you can find an appropriate mmCif file, toppology and parameters
   can be generated using eginput/PSF_generation/genLigandCif.py.

2) ACPYPE
  http://bio2byte.be/acpype/

3)  The Dundee PRODRG2 Server
  http://davapc1.bioch.dundee.ac.uk/cgi-bin/prodrg

After your PSF is read (via e.g. protocol.initStruct) or
generated (via e.g. psfGen.seqToPSF or protocol.loadPDB), use

psgGen.cisPeptide(startResid, segName)

to make a cis peptide bond between residues startResid and
startResid+1. Careful! By this definition startResid is one less than
the residue you'd like to be a cis-residue.  The segName argument can
be omitted if there is a single segid in the structure. After calling
cisPeptide, you can write out the updated PSF, which can be used in
the future without another call to cisPeptide.

  Here's a script which adds unknown atoms given some known positions.
I'l send the resulting pdb file separately.

Charles

- -- 
Charles Schwieters     email: Charles.Schwieters@nih.gov
		       www:   http://schwieters.org/cds
phone: (301) 402-4914  FAX:   (301) 402-2867

rtf @TOPPAR:topallhdg_new.pro end

evaluate ($kbbang = 500.0)
evaluate ($kbbimp = 500.0)

!
! building-on sidechains, given known bb atoms (from pdb file only).
!
! to use: change the line with input.pdb appropriately. 
! Output will be in out.pdb
!
!  set constraints to only calculate necessary interactions 
!    (don't calc bb-bb interactions)
!  don't make the improper interactions too strong: will torque bonds, angles
!
! to try: turn off vwd entirely in stage1!
!
! script based on Nilges et.al. script for random starting coords
!
! CDS 7/8/02
!


param @TOPPAR:parallhdg_new.pro end





segment
   name="    "
   SETUP=TRUE
   chain
LINK PEPP    HEAD - *     TAIL + PRO     END  { LINK to PRO }
LINK PEPT    HEAD - *     TAIL + *       END

FIRSt PROP                TAIL + PRO     END { nter for PRO }
FIRSt NTER                TAIL + *       END

LAST  CTER   HEAD - *                    END


SET ECHO=TRUE END

coor @input.pdb
end
end


delete select (name hd1 and resname his) end

!remove all hydrogens!
delete select (hydro) end

write psf output=l.psf end

param
!  @TOPPAR:parallhdg_new.pro
   bond     CP     NH3     $kbon     1.473 ! CP - proline at N terminus
   angle    CP   CP    NH3    $kang 103.2  ! CP - proline at N terminus
   angle    HA   CP    NH3    $kang 109.5  ! CP	- proline at N terminus
   angle    CT   NH3   CP     $kang 112.0  ! CP - proline at N terminus
   improper HC  HC  CT  CP  50.0 0 -70.0   ! N-term Pro
   improper  HA   HA   CP   NH3      $kchi   0   -70.874  ! CP methylene
end

coor @indira.pdb

vector do (fbeta=10) (not known)    {*Friction coefficient for MD heatbath.*}
vector do (mass=100) (not known)            {*Uniform heavy masses to speed*}
vector do (fbeta=0) (known)    {*Friction coefficient for MD heatbath.*}
vector do (mass=0)  (known)              {*Uniform heavy masses to speed*}

                                            {*molecular dynamics.          *}

constraints fix (known) end
vector do (vx = 0) (known)
vector do (vy = 0) (known)
vector do (vz = 0) (known)

evaluate ($init_t = 1000 )    {* Initial simulated annealing temperature.*}

vector do (vx = maxwell($init_t)) (attr mass>0)
vector do (vy = maxwell($init_t)) (attr mass>0)
vector do (vz = maxwell($init_t)) (attr mass>0)
vector do (x=(random()-0.5)*20) (attr mass>0)
vector do (y=(random()-0.5)*20) (attr mass>0)
vector do (z=(random()-0.5)*20) (attr mass>0)

!try bonds first
flags exclude * include bond  end

constraints 
 interaction (attr mass>0) (attr mass=0) 
 interaction (attr mass>0) (attr mass>0)
end

!mini powell
!  drop=1e5
!  nprint=1
!  tolg=1e-5
!  nstep=1000
!end

flags exclude * include bond angle dihe cdihe impr vdw  end

!energy end
!mini powell
!  drop=1e5
!  nprint=1
!  tolg=1e-5
!  nstep=1000
!end

    evaluate ($kbon = 0.00005  )                 {* Bonds.                 *}
    evaluate ($kang = 0.00005  )                 {* Angles.                *}
    evaluate ($timestep = 0.04)
    evaluate ($nstep = 100)

    !constraints 
    !	 interaction (all) (all) 
    !	 weights bond $kbon angl $kang impr $kimp vdw 0 elec 0 end 
    !end
    !dynamics  verlet
    !	   nstep=5000   timestep=$timestep   iasvel=current   
    !	   tcoupling=true tbath=$init_t  nprint=50  iprfrq=0 
    !end	    

    while ($kbon < 0.01) loop stage1
      
            evaluate ($kbon = min(0.25, $kbon * 1.25))
            evaluate ($kang = $kbon)
            evaluate ($kimp = $kbon/10)
	   
            noe scale * $knoe end
            !restraints dihed scale 0. end 
            constraints 
                interaction (attr mass>0) (attr mass=0) 
                weights bond $kbon angl $kang impr $kimp vdw 5e-4 elec 0 end 
                interaction (attr mass>0) (attr mass>0) 
                weights bond $kbon angl $kang impr $kimp vdw 5e-4 elec 0 end 
            end
	    	    
            dynamics  verlet
                  nstep=$nstep   timestep=$timestep   iasvel=current   
                  tcoupling=true tbath=$init_t  nprint=50  iprfrq=0 
            end	    
	    
    end loop stage1

    parameter                {* Parameters for the repulsive energy term.  *}
      nbonds
        repel=0.9            {* Initial value for repel - modified later.  *}
        nbxmod=-3            {* Initial value for nbxmod - modified later. *}
        wmin=0.01 
        cutnb=4.5 ctonnb=2.99 ctofnb=3. 
        tolerance=0.5 
      end
    end
    ! add vdw and slowly increase its weight
    parameter nbonds 
         atom cutnb 100 tolerance 45 repel=1.2 
         rexp=2 irexp=2 rcon=1.0 nbxmod 4
    end end
    flags exclude * include bond angle impr dihe cdihe vdw end

    constraints 
        interaction (attr mass>0) (attr mass=0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.002 elec 0 end 
        interaction (attr mass>0) (attr mass>0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.002 elec 0 end 
    end
    dynamics  verlet
          nstep=500   timestep=$timestep   iasvel=current   
          tcoupling=true tbath=$init_t  nprint=50  iprfrq=0 
    end	    
    constraints 
        interaction (attr mass>0) (attr mass=0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.005 elec 0 end 
        interaction (attr mass>0) (attr mass>0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.005 elec 0 end 
    end
    dynamics  verlet
          nstep=500   timestep=$timestep   iasvel=current   
          tcoupling=true tbath=$init_t  nprint=50  iprfrq=0 
    end	    
    constraints 
        interaction (attr mass>0) (attr mass=0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.01 elec 0 end 
        interaction (attr mass>0) (attr mass>0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.01 elec 0 end 
    end
    dynamics  verlet
          nstep=500   timestep=$timestep   iasvel=current   
          tcoupling=true tbath=$init_t  nprint=50  iprfrq=0 
    end	    
    constraints 
        interaction (attr mass>0) (attr mass=0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.5 elec 0 end 
        interaction (attr mass>0) (attr mass>0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.5 elec 0 end 
    end
    dynamics  verlet
          nstep=500   timestep=$timestep   iasvel=current   
          tcoupling=true tbath=$init_t  nprint=50  iprfrq=0 
    end	    

!flags exclude * include bond angle end

energy end
mini powell
  drop=100
  nprint=1
  tolg=1e-5
  nstep=10000
end

    constraints 
        interaction (attr mass>0) (attr mass=0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.1 elec 0 end 
        interaction (attr mass>0) (attr mass>0) 
        weights bond $kbon angl $kang impr $kimp vdw 0.1 elec 0 end 
    end

mini powell
  drop=100
  nprint=1
  tolg=1e-5
  nstep=10000
end

write coord output=out.pdb end

stop

so far I have been using standard XPLOR protocols in XPLOR-NIH. I'd like
to try the torsion angle dynamics, but have little idea of parameter
settings. Would you have an example script, such as that described in
Schwieters and Clore (2001) J. Magn. Reson., 152, 288-302, for protein
G?

Please contact Charles@Schwieters.org for that exact script. Examples of
using the IVM can be found in multiple scripts located in the eginput
directory of the Xplor-NIH distribution. In particular, see
eginput/protG/anneal.inp

These messages are normal, and should not be of concern unless only a
very small number of steps are being taken.

Using the radius of gyration restraint is similar to using an
experimental term.
First, you need to set up the restraint, and then you have to turn on
its flag.

To set up the restraint, use the COLLapse statement at the top of your
script, near where you read in the NOEs:

collapse
    assign (aSelection) 100.0 aTargetValue
    scale 1.0
end

To turn it on, include COLL in your flags statement:

flags exclude * include bond angle improper ... noe coll end

You need to keep in mind a few things in order to set up the restraint
properly:

1.  The selections should not include any flexible regions.

2.  The target value (in Angstroms) should be (2.2 * Nresidues^0.38) - 1
where Nresidues is the number of residues in the selection.

3.  For anisotropic structures, you need to divide up your structure
into several roughly spherical pieces and create a separate COLLapse
ASSIgn
statement for each.

Other than those caveats, the Rgyr restraint is relatively
bulletproof--it's
not particularly sensitive to the force constant or whatever.  And it
doesn't
affect sampling or convergence.

            Hi Dr. Lecomte,

Thanks for your email.

The "raw" DELPHIC torsions code you're trying to use has been deprecated,
in favor of using fitted Gaussians.

To use the "raw" DELPHIC torsions, you'll need to change all your
RAMA commands to XRAM. This includes every call to RAMA inside
of the setup_longrange_4D_new.tbl script.

Sorry for the incompatability--we're still cleaning up all the old code and
examples.

Note that you're almost certainly better off using the stuff in
databases/torsions_gaussians or databases/torsions_quarts, since those
approaches give *much* smoother forces on the atoms. See
J Magn Resn 146 (2): 249-254 Oct 2000 .

Please let us know if you find any more problems!

--John Kuszewski

> calculations is negative. According with the original article E is
> defined as
> - - -k*log(Pi) where Pi is the probability, that I think should be
> between 0 and
> 1, so I don't understand how E can be negative. Do you know?
> Thanks again,
> Hugo

The formulation for the DELPHIC torsion and DELPHIC position terms'
potentials of mean force is

E = -ln Pi, where Pi is the probability of being in bin i.

Pi is now defined as (nExamples(i) / volume(i)) / overallDensity.

This change was necessitated by some of the averaging I do to get a
defined potential energy
in unpopulated regions.  If there are no examples in a particular bin,
I start including that bin's
neighbors along each axis until I reach a certain minimum number of
examples.  Then I divide
by the total volume of the bins I just searched to get a "local
density."  In order for the units to
work out correctly, I need to divide this local density by the overall
density (totExamples / totVolume).

The upshot is that regions of space that have greater-than-average
density of examples now have
negative energies, and regions with less-than-average density of
examples now have positive energies.

Add a NCS-style positional symmetry term to keep the subunits
identical, and a distance symmetry term to maintain C2 symmetry.

To a Python script, one would add 

# this assumes the two subunits have segids A and B:

from distSymmTools import create_DistSymmPot, genDimerRestraints
from selectTools import minResid, maxResid

dSymm = create_DistSymmPot('dSymm',
      	                   genDimerRestraints(segids=['A','B'],
                                              resids=range(minResid(),
                                                           maxResid()+1,10)))
potList.append(dSymm)

# ``NCS'' term - keep dimer subunits identical
from posDiffPotTools import create_PosDiffPot
ncs = create_PosDiffPot("ncs","segid A","segid B")
ncs.setScale(10)
potList.append(ncs)

> 
> Could you tell me what the exact meaning of the coefficients in SANI are?
> I am confused about third coefficient which you call "rhombicity". Is this
> "R" which is "Ar/Aa", or "Ar" itself?
> 

this should answer your question: The SANI COEFficient statement takes
three arguments:

 COEFficient   

 where the three real values are DFS, anisotropy, rhombicity. Their
 meaning becomes clear in the expression for the dipolar coupling shift:

 DFS + anisotropy*( (3*cos(theta)^2-1)+
		     rhombicity*(3/2)*sin(theta)^2*cos(2*phi) )

[note carefully the parentheses.]

> We are trying to incorporate multiple dipolar coupling constraints into our
> structure calculations using xplor-nih.  We have a sample input file for
> xplor-nih that has a few force constants (CACO, NCO, HNC, CH) that appear to
> be scaled relative to the N,NH force constant (for example, $k_sani_CACO =
> 0.035*$ksani), but we are unclear about where these force constant estimates
> come from, and how to estimate the missing ones (in particular, HNCA).  Any
> help in this regard would be most appreciated.
> 

The force constants are based on the measurement error, and in the
example you ask about, the dipolar couplings are all normalized
relative to NH (which makes it easier to estimate the appropriate
force constants). 

e.g. let us say the measurement error for the NH dipolar couplings is 1 Hz,
then the final force constant should be around 1.
Let us say that the error for the NC' coupling is around 0.5 Hz.  The NC'
couplings are 8.33 times smaller than the NH ones
(from the bond lengths and gyromagnetic ratios), so when normalize to the NH
couplings (by multiplying their values by 8.33)  the error would be around 4
Hz. Since we minimize the square of the difference between observed and
calculated, this would translate into a force constant that should be 1/16
that use for the NH couplings.

question:
> I then ran the xplor and the .out, .pdb, and .psf files were
> generated. I can view my protein in the vmd-xplor (.pdb). The only
> problem is that when I try to see these files in a text editor this is
> the message that I receive:
> 
> [binary file or file with \000 character]

It's probably due to bogus characters Xplor-NIH prints when it can't
determine your username. This will be fixed in the next release. Until
then you might fix your files using the following command:

perl -i -pe 's/\000//g' file

If you think that Xplor-NIH should run on your computer, but it gives
a message like: xplor not configured for platform: $ARCH

Try adding a line to the file arch/equivList of the form
OARCH ARCH

where OARCH is the architecture string given in the distribution and
ARCH is the string which Xplor-NIH claims is not supported.

>         I have a quick question to ask you about the random number 
> generation SEED in XPLOR-NIH.  The scripts supplied in the eginputs use 
> either   $seed=778, or $seed=65748309 for the initial velocities.  Are 
> these differences related to the architecture of the computer used for 
> calculations?  I'm running Mandrake 9.0 linux, and have had to adjust 
> things like the maximum  for random number generation in other types of 
> programs (I typically use values for this of  2144563822.0), and use 
> seed values of say  82364 (which is actually what I use for CNS 
> calculations). 

The exact value of the seed shouldn't matter too much- but it
shouldn't be too large- because the XPLOR interface stores the seed as
a floating point number, strange things happen it the seed gets larger
than about 10^8.

>      Will altering the SEED change things much? Or will it 
> just produce a slightly different set of results.

The value of the seed controls the exact sequence of pseudo-random
numbers used for things like setting initial random velocities. Thus,
changing the seed will change the details of a particular
structure. However, the set of final structures should be relatively
insensitive to the seed value. If the set of results does depend
sensitively on seed, it usually means that the optimization protocol
needs some further work.

You are probably redirecting print or display output to a file. When
you do this, and you also want some output to come out stdout, you
must perform two steps: close the file, and redirect output to
OUTPUT. Here's an example:

  set print=noe.out end
  print noe thresh=0.5 end
  close noe.out end
  set print=OUTPUT end

If you do not both close the file and redirect output back to stdout,
you will run out of files if you run many structures.

Another effect of not redirecting print or display output back to
OUTPUT is that subsequent output to stdout will not occur. Rather the
data will be written to a file named fort.# (the filename is system
dependent). 

This is usually due to lines which are too long. The maximum line
length in the XPLOR interface of Xplor-NIH is 132 characters. You can
work around the error by inserting carriage returns in the offending
line, such that it is no longer too long.

This is due to a conflict with Security Enhanced Linux (SELinux) as
configured in some distributions. You can get Xplor-NIH running in one
of two ways:

   * Change the default security context for Xplor-NIH by issuing the command:

       find . -name \*.so -exec chcon -t texrel_shlib_t {} \; 

     in the top level Xplor-NIH directory.

                    or

   * Disable SELinux altogether by setting the line

       SELINUX=disabled

     in your /etc/sysconfig/selinux file. Then reboot.

Yes. By the Darwin_8_x86 packages on the download site.

Yes. In the xplor-nih distribution directory, please execute these two
commands: 

echo "Darwin_8_x86 Darwin_9_x86 Darwin_10_x86" >> arch/equivList
rm bin.Darwin_8_x86/libgcc_s.1.dylib

Applies for Xplor-NIH versions prior to 2.15:

This can be worked around by editing bin/xplor.in in the distribution
directory and removing the following lines:

if [ $arch = "Linux_2.4_i686" ]; then
    #
    # kludge for Redhat kernels with back-ported thread library.
    #
    LD_ASSUME_KERNEL="2.4.19"
    exported="$exported LD_ASSUME_KERNEL"
fi

This workaround also seems to help for newer versions of the Suse, and
Gentoo distributions.

  Scripts generating multiple structures can be conveniently
parallelized using the -parallel command-line switch. An example of
this use is shown on page 43 of the slides at 
http://bit.niddk.nih.gov/xplor-nih/nih/class08Dec2006.pdf

  Also, shared memory parallelism can be employed when performing
refinement on an ensemble of structures. Please enquire for more
details.

Please see http://bit.niddk.nih.gov/xplor-nih/parallel.txt (or
parallel.txt in the Xplor-NIH distribution directory) for details. 

An example PBS script can be found at
http://bit.niddk.nih.gov/xplor-nih/nih/xplor.pbs